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I .   INTRODUCTION 
A .   BACKGROUND 

The  Marcum-Swerling  (M-S)  models  are  the  most  commonly 
used  target  models  in  modern  radar  analysis.  The  other  target 
models  have  also  developed  over  period  of  time,  such  as  chi- 
square  (Weinstock) ,  log-normal.  However  in  this  thesis  only  M- 
S  models  will  be  discussed. 

Histrocially,  the  earliest  descriptions  of  a  target  were 
in  terms  of  a  single  cross-sectional  area  value.  This  quantity 
was  usually  some  type  of  an  average  cross  section  over  the 
aspect  angles  which  the  system  designer  considered  most 
probable.  The  approach  led  the  system  designers  to  associate 
a  single  cross-sectional  area  of  unvarying  value,  with  a 
target  and  thus  led  to  the  model  of  a  steady-state  target. 
Following  this  concept,  systems  were  designed  to  achieve  a 
somewhat  arbitrarily  specified  signal-to-noise  ratio  at  some 
specified  maximum  operating  range.  The  inherent  basic 
assumption  was  that  the  required  detection  capability  and 
parameter  estimation  accuracy  could  be  achieved  if  this  goal 
was  met.  The  dominant  rule  of  the  human  observer  in  early 
detection  systems  tended  to  make  this  approach  acceptable.  As 
the  requirements  on  radar  systems  became  more  demanding,  and 


as  automatic  processing  of  radar  returns  was  developed,  the 
situation  underwent  a  change.  The  desire  to  optimize  radar 
designs  established  the  need  for  more  precise  target  models . 

The  initial  work  of  J.I.Marcum  of  the  RAND  Corporation  in 
the  late  1940s  applied  the  work  of  Rice  of  the  Bell  Telephone 
Laboratories,  relating  to  steady-state  signals  immersed  in 
noise,  to  the  problem  of  radar  signal  detection.  Marcum 
effectively  gave  a  complete  treatment  to  the  statistical 
problem  of  a  group  of  constant  amplitude  signal  pulses  in  the 
presence  of  noise.  His  work  resulted  in  an  evaluation  of 
previously  untabulated  functions  and  a  direct  application  of 
the  results  to  a  gamut  of  problems  in  automatic  detection. 
The  statistical  approach  used  in  Marcum ' s  studies  is  the  basis 
for  much  of  the  work  that  followed  in  the  radar  detection 
area . 

The  detection  of  signals  resulting  from  a  fluctuating 
target  is  basically  different  than  for  signals  resulting  from 
a  non- fluctuating  target.  A  requirement  remained,  therefore, 
to  produce  the  same  type  of  analysis  and  supporting  tabulation 
of  basic  functions  for  a  fluctuating  target  as  had  been  done 
by  Marcum  for  the  nonf luctuating  target.  This  work  was  done  by 
Swerling,  in  the  early  1950s.  He  extended  Marcum 's  approach  by 
employing  four  target  models  and  two  different  density 
functions  in  conjunction  with  two  extremes  of  correlation. 

Swerling' s  case  1  and  case  2  are  the  target  models  which 
describe  large  complex  targets,   such  as  aircraft,   rain 


clutter,  and  terrain  clutter.  They  represent  two  extremes  of 
correlation,  and  the  statistical  model  is  derivable  from  the 
physical  scattering  characteristics  of  such  bodies.  It  can 
also  be  stated  that  empirically  derived  data  on  such  targets 
are  in  very  good  agreement  with  those  obtained  from  the 
mathematical  model . 

At  the  time  Swerling  was  doing  his  work,  it  was  realized 
that  the  model  suitable  for  large  complex  targets  would  not 
give  an  adequate  description  of  large  simple  structures .  The 
problem  of  selecting  a  model  to  describe  this  type  or  class  of 
target  is  a  difficult  one  to  handle  directly.  There  is  no 
obvious  parallel  development  for  it  from  the  scattering 
characteristics  of  a  specific  body  type.  The  approach 
employed  was  to  establish  qualitatively  the  basic  statistical 
behavior  of  the  target  cross  sections  of  interest.  Having  done 
this,  a  convenient  form  of  distribution  was  selected. 

Swerling  picked  a  class  of  distributions  know  as  the  chi- 
square  class,  which  has  the  exponential  density  function 
distribution  as  one  member.  By  an  appropriate  selection  of  a 
class  parameter,  namely,  the  number  of  degrees  of  freedom,  the 
qualitative  properties  desired  for  the  distribution  associated 
with  large  simple  structures  was  obtained.  Swerling 's  cases  3 
and  4  represent  targets  that  behave  as  if  they  have  four 
degrees  of  freedom  and  are  valid  for  targets  such  as  rockets, 
missiles,  and  space-based  satellites. 

Swerling' s  two  classes  of  density  functions  evaluated  at 


two  extremes  of  correlation,  together  with  Marcum's  constant 
target  model  (case  5),  have  been  the  bases  of  virtually  all 
radar  detection  analyses. 

Five  target  models  according  to  Marcum-Swerling  scheme  are 
as  follows : 

a)  Swerling  case  1  -  The  echo  pulses  received  from  a 
target  on  any  one  scan  are  of  constant  amplitude  throughout 
the  entire  scan  but  are  independent  (  uncorrelated  )  from  scan 
to  scan.  This  assumption  ignores  the  effect  of  the  antenna 
beam  shape  on  the  echo  amplitude.  The  probability  density 
function  of  the  radar  cross  section  a  is  given  by  the 
exponential  density  function. 

w(a,o)  =±-e(-ar3) ,  a>0  (1) 

a 

where  <5  is  the  average  radar  cross  section 

b)  Swerling  case  2  -  With  the  same  density  function  as 
case  1  but  the  f luctulations  are  more  rapid  than  in  case  1  and 
are  taken  to  be  independent  from  pulse  to  pulse  rather  than 
from  scan  to  scan. 

c)  Swerling  case  3  -  The  fluctuation  is  assumed  to  be 
independent  from  scan  to  scan,  as  in  case  1,  but  the 
probability  density  function  is  given  by  the  chi-square 
distribution  with  four  degrees  of  freedom. 

v(0'^  =  nrr\  ,  =  (^/o)*"1  e"(*^-iS  e"2^ ',  (tf=2)  (2) 

\J\-±)     !     Q  q* 


d)  Swerling  case  4  -  With  the  same  density  function  but 
the  fluctuation  is  pulse-to-pulse  according  to  Eq(2). 

e)  Marcum's  model  (case  5)  -  With  the  constant  density 
function  which  represents  steady-state  target. 

The  probability  density  function  of  equation  (1),  applies 
to  a  complex  target  consisting  of  many  independent  scatterers 
of  approximately  equal  echo  areas.  The  probability  density 
function  assumed  in  case  3  and  4  is  more  indicative  of  targets 
that  can  be  represented  as  one  large  reflector  together  with 
other  small  reflectors. 

There  are  curves  available  that  can  be  used  to  calculate 
the  probability  of  detection  for  each  of  Swerling  cases,  but 
for  parameters  between  those  charted,  the  designer  has  to 
interpolate.  This  can  be  inaccurate  sometimes.  Therefore  it  is 
convenient  and  useful  to  provide  accurate  programs  to 
calculate  for  any  parameters  needed  by  the  user. 

B.   RADAR  RECEIVER  MODEL 
1 .    DESCRIPTION 

The  typical  super  heterodyne  radar  receiver  and  the 
mathematically  receiver  model  is  depicted  in  Fig  1.  The 
difference  between  square-law  and  linear  envelope  detector  is 
that  a  square-law  envelope  detector  is  used  in  the  small 
signal  optimum  receiver  and  a  linear  envelope  detector  is  used 
in  the  large-signal  receiver.  For  mathematical  convenience, 
the  square-law  detector  is  applied,  but  the  performance  of 


both  detectors  is  practically  the  same  for  all  signal-to-noise 
ratios . 
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Figure  1:  Radar  Receiver  Model 

The  multiple-pulse  detection  is  used  to  improve  the 
detectability  of  the  signal  by  achieving  integration  gain.  The 
noncoherent  integration  provides  an  integration  gain  even  when 
the  signal  has  a  random  phase  and  is  rapidly  fluctuating,  such 
as  Swerling  case  2  and  4.  By  contrast,  coherent  integration 
needs  a  nonf luctuating  or  slowly  fluctuating  signal  with 
predictable  phase  characteristics. 


2 .   ASSUMPTIONS 

The  assumptions  embodied  in  the  M-S  detection  problems  are 
as  follows: 

•  A  received  pulse  train  consisting  of  N  samples  of  noise 
or  signal-plus-noise  is  available. 

•  The  signal  is  imbedded  in  white  Gaussian  noise  of  known 
spectral  density  of  N0/2  . 

•  The  signal  is  of  unknown  phase  type,  where  the  RF  phase 
between  pulses  in  the  train  are  randomly  distributed. 

•  The  processing  consists  of  a  matched  filter  and  square- 
law  envelope  detector  which  operates  on  each  pulse  of  the 
train  and  a  linear  integration  which  combines  the  n  square-law 
envelope  detected  pulse. 

•  The  directivity  pattern  of  the  antenna  is  rectangular  so 
that  the  n  return  pulses  resulting  during  the  antenna's  dwell 
time  on  the  target  are  unaffected  by  the  antenna's  radiation 
pattern. 

•  The  target's  cross  section  can  be  described  by  a  chi- 
squared  distribution  with  2k  degree  of  freedom. 

wio.o)    -  -j-\-   4  (4?)"1  e-i«"*  (3) 

[K-±) !   o     q 

C.   RADAR  DETECTION  PHILOSOPHY 

Radar  detection  is  complicated  by  the  fact  that  the  target 
cross  section  a  is  a  random  variable  fluctuating  with  time  and 


both  noise  and  clutter.  Extraction,  or  parameter  estimation, 
is  likewise  complicated  by  random  fluctuations  of  the  target 
echo.  For  the  radar  detection  case  Neyman- Pearson  criterion  is 
used  which  needs  neither  a  priori  probabilities  nor  cost 
estimates.  In  radar  terminology  whose  objective  is  to  maximize 
the  probability  of  detection  for  a  given  probability  of  false 
alarm.  This  objective  can  be  accomplished  by  using  a 
likelihood  ratio  test.  Specially,  there  exists  some 
nonnegative  number  r\  such  that  hypothesis  Hx  (i.e., target  is 
present)  is  chosen  when 

f     (v) 

A(y)  -    r,  \   *  n  (4) 

and  hypothesis  H0  (no  target  present)  is  chosen  otherwise. 
There  are  two  types  of  errors  which  can  be  made.  If  we  say  a 
target  is  present  when  in  fact  it  is  not,  an  error  of  the 
first  kind  is  made  (see  Fig  2)  .  That  is,  we  choose  Hx  given 
that  H0  is  true  which  is  denoted  as  the  probability  of  false 
alarm.  Similarly  an  error  of  the  second  kind  is  made  when  H0 
is  chosen  and  in  fact  Hx  is  true  which  is  probability  of  miss. 
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Figure  2:  Error  of  detection 

To  optimally  detect   the  signal  with  unknown  phase, 
consider  the  narrowband  signal  of  duration  T  given  by 


s{C)  =Aa  (t)  cos  [wQt+Q  (t)  +(J>] 
=SjCOs<|>  -s^sincj) 


(5) 


where  Sz=Aa(t)  cos  [w0t+6  (t)  ]  ,  Sq=Aa(  t)  sin[w0t+6  ( t)  ]  ,  A  is  signal 
amplitude,  a(t)  is  an  envelope  function  of  duration  T,  9(t)  is 
the  signal  phase  modulation.  Waveform  r(t)  is  assumed  to  be 
prefiltered  by  an  ideal  low-pass  filter  that  is  distortionless 
within  (-wc,wc)  and  zero  outside  the  interval;  therefore  the 
filter  will  not  affect  the  band-limited  input  signal  s(t)  but 
results  in  a  band-limited  statistical  process.  Under     this 


assumption  the  input  waveform  r(t)  can  be  described  by  the 
sampling  theorem  as 

r^t)    =  sJ(t)  +  ^(0  i  =  1,2,3,  ...  ,2fcT         (6) 

Where  2fcT  is  the  number  of  samples,  and  t=i/2fc. 
The  hypothesis  testing  can  be  described  as: 

H,  :  r  =  8  +  n        (  target  present  ) 

(7) 
HQ   :  r  =  n  {  target  absent  ) 

Then  under  each  hypothesis  the  probability  density  functions 
can  be  written  as: 


fs+n(r)= e         2a°        /  fn(r)= 1 e  2o" 


The  likelihood  ratio  A(r)  is  given  by 


(8) 


2fcT      ,  , 


««    -  ^  ■  Jfr  ^  ■   "* 


exp 


i  rT, 


exp<--±-f    [r(t)-s(t)]2  dd 
N0Jo 


(9) 


expi-4-f  r(t)2dd 
W0Jo 

=  exp[-±fT  s2(t)   dt   +  |f:r(t)s(t)  dt] 

where  E  =  J0T  s2(t)dt  is  the  energy  of  the  signal.  Substituting 
Eq(5)   into  Eq(9)  yield  the  following  expression  for  the 
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likelihood   ratio   A(r    |  <j)    ) 

A(r|<J))    =  exp{--^+-i-[rrr(t)si(t)dt]cos4) 
Nn      N„    Jo 


(10) 


-  -|[frr(t)sa(t)dt]sin(J)  } 
iVn    Jo  y 


Let    the    signal-to-noise   ratio   §    =E/N0      and 

yI[?)=-2—(Tr(t)sI(t)    dt 

y  (T)=-^L-(T  r(t)s0(t)    dt 

using  Eq(ll),  Eq(10)  can  be  written  as  follows 

A(r|(J>)    =  e"+/2  exp  i^[yx{T)   cos<t>  -  yQ{T)   sin<t>] 

=  e"*/2  exp  {yAJF  r  (T)  cos  (4>+<x) } 


(11) 


(12) 


where  a=tan_1  (yQ/yx)  ,  r(  T)  =  [yQ(  T)  +yx  (  T)  ] 1/2  is  the  envelope  of  the 
radar  signal  out  of  a  filter  matched  to  the  waveform  (2/N0 
\\f1/2 ) s  (t)  ,  sampled  at  time  t.  Assume  ({)  is  a  uniformly 
distributed  density  function  U(0  2k) .  Averaging  with  respect 
to  <j)  yields  the  averaging  likelihood  ratio  as  follows 

A(r)  =e*/2—  f2nev^r(Dcos(4.*o)d<j) 

271  JO 

(13) 

=e-*/2J0[vAJFr(r)] 

where  I0  (\\f1/2r  (T) )    is  the  modified  Bessel  function  of  the  first 
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kind  and  order  zero.  Therefore  the  threshold  test  becomes 


H, 


J0[v^  r(fl)]  >  e*'2  A(r) 


(14) 


Hn 


where  ey/2A(r)  is  the  operating  threshold.  The  corresponding 
receiver  implementation  is  shown  in  Fig  3 
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Figure   3 :    Optimum  detector  for  radar  signal  of  unknown  phase 

D.      OBJECTIVE 

The  Mar  cum- Swer  ling  target  model  is  the  most  commonly 
used  model  to  calculate  radar  range  performance.  The  objective 
of  this  thesis  is  to  review  the  underlying  theory  of  radar 
detection  for  this  model  and  then  develop  a  MATLAB  programs  to 
compute  probability  of  detection  and  maximum  detection  range. 
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II.  MODELS  AND  MATHEMATICAL  METHODS 

The  procedure  for  the  calculation  of  average  probability 
of  detection  PD  for  Marcum  and  Swerling  cases  is  as  follows: 
DFind  the  single  pulse  characteristic  function,  C1(p)  ,  by 
transforming  the  single  pulse  probability  density  function 

fs+N(y)  • 

2)Find  the  impulse  characteristic  function,  CN(p)  =  [C1  (p)  ]N . 
3) Average  the  impulse  characteristic  function  over  the  target 
distribution  to  find   the  average  of  CN(p)  (    CN(p)  )  . 
4) Transform  CN(p)    to  find  the   average  ensemble  probability 
density  function  fs+n(y)    and  fn    (y)    for  the  N   pulse  return. 

5)  Find  probability  of  false  alarm  P£a  by  integrating  fn(y) 
from  yb   to  infinity. 

6)  If  Pfa  is  given  find  yb  by  means  of  the  mathematical 
recursive  method. 

7) Find  the  average  probability  of  detection  Pd,  by  integrating 
the  density  function  of  signal  plus  noise  fs+N(y)  from  the 
threshold  yb   to  °o  . 

A.   MARCUM' S  (NON- FLUCTUATING)  STEADY-STATE  TARGET  MODEL 

The  non-fluctuating  target  model  is  applied  to  spherical 
or  nearly  spherical  objects,  such  as  balloons,  many 
wavelengths  in  diameter.  The  target  model,  therefore,  can  be 
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represented  by  a  constant-valued  radar  cross  section. 

The  Rayleigh  and  Rician  probability  density  function  of 
the  envelope  detected  output  for  N  pulses  are  given  by 


fjr^)    =  rie-ri   ;    r^O 


-ul+t) 


(15) 


*B*nl*iW    =  rie       2   V^iV^)  ;  r^O 


where  \\i=E/N0  is  the  single  pulse  peak  signal-to-noise  ratio 
for  a  steady  target.  Since  the  noise  received  in  the  ith 
observation  interval  is  assumed  to  be  statistically 
independent  of  the  noise  in  every  other  observation  interval, 
and  the  signals  are  noncoherent  between  the  different 
observation  interval,  the  initial  phase  (f^  in  Eq(5)  in  the  ith 
interval  is  also  statistically  independent  of  <f>j  for  all  i^j  . 
The  likelihood  ratio  test  in  Eq  (9)  can  be  written  as: 

A(r)-fliri'r2'  I  '  •z»\8x,s2...sN)    m    fi  £iUijfi) 
foUTfO)  i-i  fotrJO) 

=  ne-*/2J0(riV^r)  <16> 

i-l 

=  e-^/2n  J0(riV^JF) 

i-l 

From  Eq(16),  the  test  statistics  is 


rinJ0(riV^r)  >  e"2  A(r)  (17) 
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For  small  signal,  the  modified  Bessel  function  I0(x)    can  be 
approximated  for  x<l    as: 


I  (x)    =  l  +  i£i+^.+ 
0  4   64 


(18) 


X2   ,  X4    .     x  _  X2   .  ^/  __4. 


lnJQ(x)  -  in(l  +  — +^-  +  .  .  .)  *^_+f  (x 

0  4   64         4 


Therefore  Eq(17)  can  be  rearranged  and  modified  as 


N 


£rht,'  d9) 


where  r|'  is  the  operating  threshold  which  is  determined  by 
specifying  the  false  alarm  probability. 

Marcum  utilizes  a  square  law  detector  law  which  allows  the 
signal-plus-noise  probability  density  function  to  be  expressed 
directly  in  the  signal-to-noise  ratio  \\f.  To  simplify  the 
calculation,  let 

N  2      N 

*-i:-T-2>  (20) 

x-i  ^    ITi 

which  Y   is  compared  to  a  suitably  modified  threshold  Yh    . 

To  find  the  probability  density  of  Y,  first  the 
probability  density  of  f5+n(qi)  is  found  by  Jacobain 
transformation  fs+n(ri)    from  Eq(15)  into  fs+n(qi)    as  follows: 
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Since  random  variable  Y  is  the  sum  of  N  statistically 
independent  random  variables  qL,  by  applying  the  relationship 
between  the  sum  random  variable  and  the  individual  random 
variable  of  the  characteristic  function,  and  the  Campbell  and 
Foster  tables  of  Fourier  transforms  ;  the  CY(p)    is  : 

cY(P)  =  n  c    (P)  =  n  r  ejpg*fs+t,(qi)  dq± 

i-1   J      i-iJ-» 

=  fi  ["  ejPQi  e~l«**/2)  IQ(JZq^)   dq,  (22) 

i-0  Jo 

-  e     e  2  dtp) 

From  Campbell  and  Foster  Tables  of  Fourier  transforms  pair 
650.0,  the  probability  density  function  of  fs+n(Y)  can  be  found 
by  taking  the  inverse  transform  of  CY(p) . 


Wfl  =  ^)W"1,/2  e"™*/2  Vi(v^AJJT),     7*0    (23) 

For  small  Y,      the  modified  bessel  function  I^fY)     can  be 
approached  as : 

j(y)=_Zf!_  [i  + £! + Z! +  ...]      (24) 

2mm\  22{m+l)      25(/t?+1)  (777+2) 


Therefore   the   resulting   normalized   square-law   detected 
probability  density  function  is 


yN-le-Y-tHf/2 

(N-l)   ! 


Wy)  =  ^;f//7  ,   y,0  (25) 
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for   noise   only    the   probability   of   density    function    f„(y)    is 


*«(«    " 


YN-le-Y 


Y*0 


(N-l)  !  ' 

A    computer    simulation    of    fn    (y)  at    the    square    law 
output    is   given    in   Figure    4 : 


(26) 
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Figure  4 :  Simulation  of  noise  at  suare  law  detector  output 

After  square-law  detection,  the  normalized  square-law 
detected  variate  (y)  is  summed  over  the  N-pulses  in  the  pulse 
train  and  then  compared  against  a  normalized  threshold  voltage 
(yb)  to  determine  the  presence  or  absence  of  a  target  return. 

The  probability  of  false  alarm  can  be  obtained  by 
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integrating  from  the  threshold  (yb)    to  infinity 


yN-le-y 


P      =   (     Y      %1  dY  (27) 

fa  Jyb    (N-l,     ' 


The  incomplete  Pearson  gamma  function  which  is  very  useful 
for  describing  M-S  model  is  defined  as: 

I(u,p)    =f^u<"+1))  ^Ilfdy  (28) 

Jo  pi 

And  in  term  of  the  Pearson's  incomplete  gamma  function,  Eq(27) 

is 

Pfa(N,Yb)    =  l-I[Yb/N^2,  (N-l)]  (29) 

where  u=(yb/N1/2)    and  P  =  N-l. 

Eq(29)  can  be  successively  integrated  by  parts  and  given 
as : 

Pfa(N,Yb)  =  if1*''* [1+JE3L+  (^-D  (^)  ,, , ,] 

*  (30) 

"-1  y/  e~r> 


K'O 

The  above  equation  also  can  be  approximated  for  N>>1  by 
letting 

Yh  yb2  N-l  (31) 

and  by   applying   Stirling ' s   approximation 

m=JZnN{2)N  (32) 

e 

The  probability  of  false  alarm  in  Eq(30)  can  be  approximated 
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as 

I — xt-   y    AT     -  Yb*N 

Pf  .OT(2j>)   e  ^  (33) 

The  probability  of  detection  for  steady-target  (case  5) 
can  be  obtain  by  integrating  fs+n(Y)  from  yb  to  °°  and  is  given 
by 

(34) 

The   Q    function   is   defined  as 

CW(«,P)=/~  v(-J)*"1  exp[-(-5!|Z!)]    i^Uvldv  (35) 

For  i\7  pulses  the  probability  of  detection  PD  can  be 
written  in  terms  of  Q  function  as  follows: 

PD(y)=Q(y/ZN$,y/ZYD  <36> 

Another  approximation  can  be  made  by  using  Gram-Charlier 
series  expansion  (  Appendix  A  )  and  noting  that  Gaussian 
family  is  closed  under  linear  operation. 


dCY(p) 
dp 


*i   =  (-J)   V       Ip-o  =  N(l+y/2) 


™2    =     (-J-)2d2^(P)  |     =  N2{1+y/2)*+N(l+y)  (3?) 

dpN 
o2=m2-ml   =  tf(l+i|r) 

where  %  is  the  ith  moment  of  the  distribution  of  Y   and  a  is 
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the  variance  of  Y. 

i      p-(r-w(i+t/2)] 


fe+n{Y)%2*N(l^)        2HT(1+*)  (38) 

f  (7)  g   1   e(-r-a)»/2tf 


Therefore  the  probability  of  false  alarm  can  be  represented 
as  : 


12M1  "   ,    z£.  Yh-N 

Pfa   «  I — —  e   2W  d7  *   /  ^^e  2  dz   where  Z=  b 


[     1     c      2N    dY  *      [       ^ 

i^^N  j^v/^  V^    (39) 


V7J 


-♦t^T] 


v/3V 
Y^  can  be   solved   from   the   above   equation 

7j,  *  v^fMPfJ  +#  (40) 

The  probability   of   detection   can  be   obtained  as 

(y^-^d+T/2))  .  .... 

Pd  s  4>    £-=== <41> 

B.   SWERLING'S  (FLUCTUATING)  TARGET  MODELS 

1.   SWERLING  CASES  1  AND  2 

The  radar  cross  section  of  a  target  (a)  is  the  area  of  a 
hypothetical  reflector  that  scatters  all  radar  beam  energy  it 
intercepts  omnidirectionally,  such  that  it  produces  an  echo  at 
the  radar  antenna  equal  to  that  energy  actually  received  from 
the  target;  that  is 
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(power  received  by  antenna) 

{incident  power  density  per  An  rad) 

(42) 


=  lim  4nR2- 


Er 


l^l2 
Where  E2  is  the  electric  field  incident  on  the  target,  ER    is 
the  reflected  field  as  measured  at  the  radar  receiver  antenna, 
and  R    is  the  distance  from  the  target  to  the  radar  receiver 
antenna . 

The  total  electric  field  received  from  a  complex  target 
can  be  expressed  as  a  summation,  namely, 


**=l  El** J  expC-^p-^) 


J47td*)  I  (43) 


and  therefore  the  radar  cross  section  is 

j4TtdK 


o  =  IE/T^T"exP(^T^)l2 
=  IE/T5^r  [cos(^-^)    +  jsin(— ,-*)]  |2 


(44) 


where 

a*  is  the  crossing  section  of  individual  scattering 
elements, 

n   the  number  of  scattering  elements, 

dK   the  range  of  the  Kth  element, and 

X    the  wavelength. 

In  order  to  relate  the  proceeding  to  a  real  target  two 
assumptions  shall  be  made: 
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(1)  On  a  short  term  basis  4KdK/X  is  a  random 
variable  which  can  take  on  any  value  from  0  to  2n  with  equal 
probability. 

(2)  The  individual  scattering  elements  have  equal 
radar  cross  sections,  that  is,  ak=c0 . 

Therefore  the  random  variables  can  be  expressed  as  the  x 
component  (ct0)  (1/2) cos  (AKdk/\ )  and  the  y  component  (o0)1/2 
cos(4Kdk/\),  the  problem  of  determining  the  probability 
density  function  of  a  is  identical  to  the  problem  of 
determining  the  distance  moved  from  an  origin  in  a  two 
dimension  walk  problem  of  n  steps,  where  the  length  of  each 
step  is  (G0)1/2  an<3  the  direction  of  each  step  is  perfectly 
random . 

The  probability  of  having  a  component  (x,y)  after  n 
steps  where  n   is  a  comparativity  large  number  is  given  by 

w{x.y)   dx  dy  -    e*p[-(*^)/*q0]  ^  dy  (45) 

nnaQ 

Converting  from  rectangular  coordinates  to  polar  coordinates 
yields 

w(R,Q)    dR  dQ  =  -JZ—exip[-(R2)  /naQ]    dR  dQ  (46) 

nna0 

The  marginal  distribution  of  R,    obtained  by  integrating  the 
preceeding  with  respect  to  9,  is 

w(p)    dR    =  —exv[-(R2)/naQ]    dR  (47) 

naQ 
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Now, since   G=R2=x2+y2  and  dG=2R   dR   Eq.(47)    becomes 

exp[-o/.no0] 


w(o)    =   — : °—  for    o  >  0  (48) 


na0 


w(G)=0  for  G<0 .  This  demonstrates  that  the  probability 
density  function  of  the  cross  section  is  exponentially 
distributed  with  an  average  cross  section  3=nG0.  Thus,  the 
average  total  echoing  area  is  the  sum  of  the  individual 
echoing  areas  of  the  individual  elements . 

a.  Case    1    Scan-to-scan    fluctuations    of    exponential 
density  function. 

The  scan-to-scan  fluctuations  of  an  exponential  nature  can 
be  applied  to  targets  such  as  jet  aircraft  when  a  radar  having 
a  fairly  high  pulse  repetition  rate  and  scan  rate  is  employed. 

For  case  1,  Swerling  derived  the  probability  of  detection 
by  initially  defining  a  target  model  where  the  received  signal 
power  is  exponentially  distributed,  namely, 

w(x,x)    =    exp  ^Lx/^  for    x  >  0 

x  (49) 

=  0  for    x  <  0 

where  x  =  input  signal-to-noise  ratio. 

x  =   average  of  x  over  all  targets  fluctuations. 
The  characteristic  equation  for  the  probability  density 
function  resulting  from  the  integration  of  N   pulse  returns 
from  an  exponential  fluctuating  target,  when  there  is  complete 
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correlation  between  pulses,  is  derived  from  the  non- 
fluctuating  case  (Swerling  case  5)  characteristic  function  as 
follows : 


iJSL) 


(p+DN 

_e 
(p+1) 


Nx(-2-) 


_       --         «    "FT  (50) 

C„(p)  -    v(x,x)  .2. —  dx 

Jo  (d+i)n 


i 


(p+1)""1  [l+p{l+Nx)  ] 

For  noise  only,  the  characteristic  equation  is  the  same  as 
for  case  5 . 

cn= ^-m  <51) 

(1+P)N 

The  probability  density  function  fn  (y)  obtained  by  using 
Campbell  and  Foster  431,  p. 44,  is  the  same  as  the  probability 
density  function  of  the  nonf luctuating  case.  i,e, 

f  (Y)  =Z f_L  (52) 

nK    '       (1+JV)  ! 

For  signal  plus  noise,  the  probability   density  function 
fs+n(Y)    was  obtained  by  using  Campbell  and  Foster  pair  581.7, 
p.  64.  For  N>1,    the  density  function  fs+n(Y)    is 
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f e+n  (Y)=[~ e j2npy  dp 


=  -L   e-^(1+AeF)    ^__ [r(tf-l)-r(tf-l,-fl&-y] 

-      „  r(jr-lf-^jr)  (53) 

wx  ak  r(w-i) 

=  _L    e-y/d*Aff)     (i+^)W"2J( £ ,JV-2) 

Wx 


where  V(v,z)  is  the  incomplete  gamma  function  and  Ifu,pjis  the 
incomplete  Pearson  gamma  function. 

For  N=l,    the  probability  density  function  can  be  obtained 
from  Campbell  and  Foster  pair  43  8,  p. 45  as 

fB.n(Y)=e^Y^/(l+x)  (54) 

therefore  the  probability  of  detection  can  be  obtained  by 
integrating  Eq(54)  from  0  to  Yb  and  the  result  is  given  by 
Swerling  as: 

l-PD=f*b  fa.a(y)    dv=I[  — p^  ,N-2]-(l+Nx)f(Yb)  (55) 

From  above,  PD can  be  obtained  as 


Pn=l-X[-^=/  (tf-2)]+(l+-i.)J" 

j( 5 ,^-2)e-V(1+>fi?) 

[l+d/i^v^FT] 


(56) 
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Another  approximation  can  be  made  by  inspecting  Eq(56)  and 
let  iVx>>l  and  Pfa<<l  such  that  the  Pearson  incomplete  gamma 
functions  are  close  to  unity.  Then  Eq(56)  can  be  rewritten 
approximately  as 

P^dt-A-J^e^fe    fl\  (57) 

Nx/2  Ffa<1 

Taking  the  logarithm  of  Eq(57)  and  using  a  series  expansion  of 
each  term, 

lnPd*(tf-l)ln(l+-4j  " — 

Nx        (Nx)   (1+-4:) 
Nx 

*(N-1)  [-r±—± 1 -1 1 -.  .  .]      (58) 

(Nx)      2  (Nx/2)2     3  (Nx/2)3 

b    Nx/2      (Nx/2)2      {Nx/2)3 
Taking  only  the  first  terms  and  using  Eq(40)  yields 

lnpD«— ^(Yt-N+l) 

Nx 

.—L[y/n$-Hp  f.)+i] 

Wx 


(59) 


Jb.  Case  2  Pulse-to-pulse  fluctuations  of  exponential 
density  function. 

Pulse-to-pulse  fluctuations  of  an  exponential  nature  apply 
to  targets  such  as  propeller-driven  aircraft  (if  the 
propellers  contribute  a  significant  portions  of  the  echo 
area) ,  to  targets  where  small  changes  in  orientation  would 
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establish  significant  changes  in  echoing  area  (such  as  long 
thin  subjected  to  a  high-frequency  signal),  or  to  targets 
viewed  by  radars  with  sufficiently  low  repetition  rates. 

For  case  two,  the  received  signal  return  still  belongs  to 
exponential  distribution.  For  this  case  the  signal  is 
completely  decorrelated  (pulse  to  pulse  fluctuations  exist) , 
the  probability  density  function  for  noise  only  (signal-to- 
noise  ratio  equals  to  0)  still  the  same  with  Case  1.  The 
probability  of  density  function  of  signal  plus  noise  for  case 
2  can  be  obtained  with  the  same  procedures,  but  have  the  pulse 
to  pulse  integration.  The  characteristic  equation  for  the 
single  pulse  is  obtained  by  letting  the  characteristic 
equation  of  case  1  to  be  N=l 

q(p)  = 1   ^  (60) 

[l+p(l+x)] 
For  impulses  which  is  completely  decorrelated  (independent), 
the  characteristic  equation  is  just  the  Nth  power  of  single 
pulse 


5"(p)  =  [^)]N -  [T^b*]"  (61) 


+p  ( 1  +x) 

Therefore  the  corresponding  probability  of  density 
function  obtained  by  using  Campbell  and  Foster  "Tajbies  of 
Fourier    transforms"   pair  431,  p. 44  yield: 
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(N-l)  ~(-y/(l+3r)) 


[ ]Ne2«pydy  =   y % (62) 

-   l+p(l+x)  (1+5^)^(^-1)  ! 

The  probability  of  detection  can  be  expressed  as: 
PD  =  i-f~fs+n(y)  dy 

(N-l)  \  Jo     1+x  1+x  (63) 


-  1-J[ % ,  {N-l)] 

(l+x)v/ft 


Another  approximation  can  be  made  by  applying  that  the 
Gaussian  family  is  closed  under  linear  operation  and  by  using 
Gram-Charlier  series  expansion  (Appendix  A)  so  that 


dCY(p) 
dp 


"k  -  (-j)     *,     U  =  ^d+^) . 


,  ...  d2CY{p)  .  .         .   .      _  ,  (64) 

^  "  (-J)2   /„   p-o  =  N[N+i)  U+3F)2. 

Variance   :  o2  =  /^-/n^/Vtl+x)  2 
The  approximate  probability  density  function  then  is 

-[Y-N{l*x))2 
f         (Y)= 1 o   2tf(l+x)2 

v^2li^(l+x)  (65) 

Therefore  the  probability  of  detection  and  false  alarm  can 
be  obtained  by  integrating  from  the  threshold  level  yb  to  °° 
and  are  represented  as : 
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..  {Yb-N(l+x)  ).  n       ..    (Yb-N)  .  ,__. 

Pd=4>[ — *=- — - — ]   ;  Pfa=b\. — 7^-]       (66) 


2.   SWERLING'S  CASE  3  AND  4 

For  Swerling's  case  3  and  case  4,  the  radar  cross  section 
can  be  described  by  chi-square  distributions  with  four  degrees 
of  freedom.  The  density  function  is  commonly  associated  with 
tabilized  missile  tankage  and  can  be  expressed  as: 

-(  — ) 
w(x,x)    =  —   e      *  x>0  (67) 

3? 

a.  Case  3  Scan-to-scan  fluctuations  with  a  chi-square 
density  function  with  four  degrees  of  freedom. 

The  characteristic  equation  for  case  3  which  represents 
the  condition  of  complete  correlation  (scan-to-scan,  no  pulse- 
to-pulse  fluctuation)  is  given  by  : 


ciP)  =  r  m  e-^  £^L  * .      ip^2z  ,         (68) 


The  characteristic  equation  and  the  probability  of  density 
function  for  noise  only  remain  the  same  as  in  case  1  and  2 

(69) 

i e^pydp   =  Z- ^_ 

-  (p+l)N  (N-l)   ! 

For  signal  plus  noise  the  probability  of  density  function 
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can  be   expressed  as 


Wy>   ■    f     (P     }  -     a  e*""cfe 

"    [1+P(1+^)1 

2  (70) 


rW-l      ^-K 


y        e  F  [2   N  y 

(N-l)l     [(l+(Nx/2)]2      '  '    1*{2/Nx) 


where  the  density  function  comes  from  Fourier  transform  pair 
581.1  in  the  Campbell  and  Foster  tables  given  as 


yo-l  G-(pa)  Y/2  x 


(s+p)a+v(s+a)a"v  r(2a)(p-a)a  (71) 

«v,a-i/2[(p-0)i1  ,Y<0 


where 


H+!   -J 


W,  (Z)  =Z  2  e  2  F^-v*  "  2n+l;Z) 


(72) 


Eq(72)  can  be  simplified  using  the  relationships  for  the 
confluent  hypergeometric  function,  that  is 

F1{2,N;Z)    =    {Z+2-N)F1{l,N;Z)+N-l 

7  (73) 

FAl.NiZ)    =  ez  Z'^iN-l)   !  I[         Z    ,N-2] 

Swerling  uses  two  identities  relating  to  confluent 
hypergeometric  functions  in  order  to  expand  the  proceeding 
into  more  familiar  Pearson's  form  of  the  incomplete  gamma 
function,  I(u,p)  which  is  defined  in  Eq  (30),  thus  Eq  (73) 
becomes 
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f    {  }  =_[iM2/^)r2z  r[ y  fN-2]  e(-y/i*(iR/2)) 

Bny         [i+(Nx/2))2  [i+2/Urx)]>/2FT 

-  (^-2)  n+ij/^l""1  jr Y. lN-2]  (74) 

[l+(Nx/2)]2  [l+(2/NJ?)]yfFFl 

e-y/n+m+  y^e-y 

[l+(Nx/2)]  {N-l)   ! 


The  probability  of  detection  can  be  obtained  by- 
integrating  the  density  function  from  the  threshold  Yb  to  °° 
as  : 

vN-2      -Yb      ov         N-2  vm  _  ~2Yb 

D     (N-2)   !  Nx+2    h  rnl  Nx 

2Yh     m  (75) 

-km  o\        2Y  N~2     ^I±   ( ^) 

{1_2(N-2)  ^_f±b_]  {1.y  e  2^      2_+NX      )]  .  N±2 

NX         Nx+2  j£o         ™* 

b.  Case  4  pulse-to-pulse  fluctuations  with  a  chi- 
square  density  function  with  four  degrees  of 
freedom 

With  no  correlation,  the  characteristic  equation  for  a 
single  hit  is  obtained  by  letting  N=l   in  Eq(52) 

C(d)    =      1+p 

r    ,   x\i2  <76> 

[l+p(l  +  |)] 

The  characteristic  equation  for  the  sum  of  N  pulse  can  be 
obtained  by  the  Nth  power  of  a  single  hit 
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CuW^C^p)]"  =   {1+P)N     2N  (77) 


l+p(l+x/2)] 

The  inverse  Fourier  transform  of  above  equation  is 
obtained  from  Campbell  and  Foster  transform  Pair  581.1, p. 64, 
this  yields  the  following  expression  for  the  probability  of 
density  function  in  the  condition  of  signal  plus  noise 

t„w  =r|e-(f — ^&N  2N 

j^3?  [l+p(l+x/2)]2N 


*-i  e  i+x/2 


y"  x  e 


(N-l)  !    F1(-N,N;      *L2  y) 


(1+X)2tl  1  l+x/2 

The   confluent  hypergeometric   function  F1(-N,N;a)   can  be 
expanded  as : 


1  N    1!     N(N+1)         2!     N{N+1)  (N+2)  3! 


(79) 


K=0 


(-1)K[N\/(N-K)   !]    sA 
[  (N+K-l)  1/  (N-l)   !]  K\ 


From  Eq(78)  and  Eq(79),  the  density  function  fs+n(Y)  can  be 
rewritten  as 


(80) 


_  y*t-ie    i»*/2#j  _   x/2  v  K yK 

S"  "   (l+x/2)2"  M   l+x/2   [(#+*-!)  !  (N- 


K)   !  JC!  ] 


The  probability  of  detection  is  obtained  by  integrating  Eq(80) 
from  the  desired  threshold  Yb  to  °°.  From  the  definition  of 
incomplete  gamma  function,  the  probability  of  detection  can  be 
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written  as 


pd  -  f"  Wy)  ^y 

(81) 


r* 


*  _     J[ — _,   ; ,iy+*-i] 

x_     JY!     yv  (XJ     (H-X/2)  y/i\F7C 

(l+x/2)*jfe*   2J         JC!  (JNT-JO  ! 


With  Eq(40)/  the  probability  of  detection  can  be  approximated 
as  : 


.  y/N   4>_1(Pfa)  +N 


P     -  i-    N[  V  {  1)K         d+q/2)  y/^X 

°      (1+0/2)  *f=j   2J         K\     (N-K)   ! 


C.  SEARCH  RADAR   DETECTION     RANGE   CALCULATION 

The  Marcum-Swerling  theory  represented  by  extensive  sets 
of  curves  from  the  computer  programs  can  be  used  to  determine 
the  detection  range  of  a  practical  radar  by  introducing 
detection  loss  and  others  parameters.  There  are  many  types  of 
detection  losses  which  have  been  identified  so  far,  and  when 
these  are  considered,  reasonable  predictions  of  radar 
performance  can  be  obtained. 

The  radar's  detection  range  can  be  determined  by  applying 
the  desired  signal-to-noise  ratio  determined  from  the  Marcum- 
Swerling  theory  and  the  calculated  detection  loss,  as  given  by 
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Pt.Gt.G-K   o    ,  ,  /.    .  . 

^ax=[ *-*-= ^-11/4   (m)  (83) 

(4k)  2kTBFnL(±) 
N 


where 

i^x  =  The  maximum  detection  range  in  meters 

Pt  =  Peak  transmitter  power  in  Watts 

A,  =  Wavelength  in  meters 

Gtr  =  Transmitter  and  receiver  antenna  power  gains 

a  =  Average  radar  cross  section  in  square   meters 

k   =  Boltzmann's  constant=  1.38  x  10"23  J/deg. 

T  =  Effective  system  input  noise  temperature  in  degrees 
Kelvin  (°  K) 

B   =  Receiver  bandwidth  in  Hertz . 

L  =  Detection  system  power  loss  factor 

Fn=  The  receiver  noise  figure 

S/N=   The  smallest  output  signal-to-noise  ratio 

When  n    pulse  are  integrated  previous  equation  can  be 
written  as 

PtGtGrK2an        .... 

^nax=[ L-^£ ^]1/4   («)  (84) 

(4n)2kTBFnL(±)n 
N 

where  the  parameters  are  the  same  as  that  of  Eq  (83)  except 
that  (S/N)n  is  the  signal-to-noise  ratio  of  one  of  the  n 
pulses  that  are  integrated  to  produce  the  required  probability 
of  detection  for  a  specified  probability  of  false  alarm. 


34 


III.   SOFTWARE  DEVELOPMENT 

This  chapter  describes  the  development  of  MATLAB  programs 
for  the  efficient  and  accurate  computation  of  probability  of 
detection  based  on  Marcum  and  Swerling  theory  of  radar 
detection.  The  MATLAB  source  code  is  given  in  Appendix  B,  and 
complete  programs  are  available  from  Professor  G.S.  Gill,  Code 
EC/G1,  Naval  Postgraduate  School,  Monterey,  CA  93943. 

A.   PROGRAM  STRUCTURE 

The  overall  program  structure  is  shown  in  Figure  (5) .  The 
structure  is  that  of  a  main  menu  program  which  calls  various 
submenu  programs  (mscurve.m  number. m  number. m  )  as  required. 
The  submenu  programs,  when  called,  will  then  display  the 
purpose  of  the  subprograms.  When  called,  the  subprograms  will 
call  the  function  programs  to  do  the  actual  computation.  It  is 
possible  to  exit  the  process  from  either  the  main  program  or 
from  the  subprograms  or  the  function  programs.  The  advantage 
of  this  format  is  that  if  the  user  wants  to  change  one  of  the 
subprograms  or  function  programs  and  wants  to  add  other 
programs,  this  system  can  accommodate  it.  For  each  subprogram 
there  is  an  error  detect  prevention  and  data  entry  double 
check  function  to  inform  the  user  and  restart  the  process. 
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[MAIN  MENU] 
MS.M 

I 

[  SUBMENU.M  ] 
MSCURVE.M 

[  SUBMENU] 
NUMBER.M 

[  SUBMENU] 
RMENU.M 

1 

IE 

I 

[ SUBPROGRAM  ] 

RADAR. M 

RADAR1.M 

MEYER.M 

[SUBPROGRAM] 

POINT.M 
RADTHRE.M 

[  SUBPROGRAM  ] 
RANGE1.M 
RANGE3.M 
RANGE4.M 
DETECT.M 
DETFRACT.M 

[  FUNCTION  PROGRAMS  ] 

SWERL1.M 
SWERL2.M 
SWERL3.M 
SWERL4.M 
SWERL5.M 

I 

[  FUNCTION  PROGRAMS  ] 

THRESH.M 
THRESHM.M 

PROB.M 

MARCUM.M 

Figure  5 :  Program  structure 

1.   User's  guide  and  instruction 

To  use  these  programs,  the  following  MATLAB  files  have 
to  be  copied  to  the  user  subdirectory.  A  brief  explanation  of 
each  is  also  given. 

•ms.m  is  the  main  menu  program  and  gives  brief  descriptions  of 
the  M-S  model  and  displays  the  main  menu. 

•mscurve.m  ,rmenu.m,  and  number. m  are  the  submenu  programs. 
These  give  the  purpose  of  the  programs  when  called. 
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•  radar. m,   radral.m,   meyer.m,   are  the  subprogram  used  to 
integrate  swerll.m,   swerl2.m,   swrl3.m,   swerl4.m,   swerl5.m 
function  programs  to  calculate  and  plot  the  curves  pd  vs  S/N, 
pd  vs  pfa  and  pd  vs  N,  respectively. 

•rangel.m,  range3 .m,  range4.m,  are  the  subprograms  responsible 
for  integrating  and  calculating  the  detection  range.  The  user 
inputs  the  specified  detection  loss  from  detect. m.  The 
rangel.m,  range3.m,  range4.m,  will  load  the  data  from  detect. m 
and  calculate  automatically  and  display  the  detection  curves 
on  the  screen. 

•point. m  radthre.m,  are  the  subprograms  responsible  for 
calculating  and  displaying  the  numerical  results  from  the 
function  program  radpoint.m  and  marcum.m,  respectively, 
•swerll,  swerl2,  swerl3 ,  swerl4,  swerl5,  prob.m,  thresh. m, 
threshm.m,  bound. m,  radpoint.m,  noise. m,  signal. m,  are  all 
function  programs  responsible  for  calculating  the  data  from 
the  subprogram  and  then  return  the  numerical  value. 

A  3  86  or  486  personal  computer  is  suggested  for  greater 
speed.  To  start  the  program  type  ms  and  press  [enter]  at  the 
MATLAB  prompt . 

2.   Printing  graphical  outputs 

Graphics  output  from  all  programs  will  be  stored  as  meta 
files  automatically.  The  operator  can  print  out  the  desired 
graphics  from  MATLAB  by  typing  Igpp  <filename>  [enter]  .  Once 
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the  user  restarts  the  programs,  the  previous  meta  files  will 
be  deleted  automatically. 

3.   Programs  options 

Once  the  ms.m  command  is  given,  the  first  screen  seen  by 
the  user  is  as  shown  in  Figure  6 . 


%THIS  PROGRAM  IS  DESIGNED  TO  ALLOW  THE  STUDENT  TO 
%VARY  THE  PARAMETERS  OF  THE  VARIOUS  SWERLING  CASES 
%IN  ORDER  TO  STUDY  THE  EFFECTS. 

CASE  #:     DESCRIPTION 

1.  Returned  pulses  are  of  a  constant  amplitude 
over  one  scan,  but  are  uncorrelated  from 
scan  to  scan. 

2.  Returned  pulses  are  uncorrelated  from  pulse 
to  pulse  and  correlated   from  scan  to  scan. 

3.  Returned  pulses  are  of  a  constant  amplitude 
over  one  scan,  but  are  uncorrelated  from  scan 
to  scan. 

4.  Returned  pulses  are  uncorrelated  from  pulse  to 
pulse  and  correlated  from  scan  to  scan. 

5.  The  static  case  with  constant  S/N  and  pulse 
amplitude 


Figure  6:  Main  menu  descriptions 

This  screen  gives  a  brief  description  of  the  five  target 
models.  After  pressing  the  "enter"  key,  the  user  will  see  a 
second  screen  as  shown  at  Figure  7 . 
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MAIN  MENU  

1)  THE  M-S  CURVES 

2)  NUMERICAL  DETECTION  PROBABILITY  CALCULATION 

3)  RANGE  DETECTION  CURVES 


Select  a  menu  number: 


Figure  7 :  Main  menu 

At  this  point  the  operator  has  different  choices  to  make 
depending  upon  what  he/she  wants  to  do.  Once  the  operator 
chooses  one  of  the  items,  the  main  menu  program  will  transfer 
to  the  selected  submenu  in  Figure  8. 


SUB  MENU  —  [THE  M-S  CURVES 

ANALYSIS]  

1)  PROBABILITY  OF  DETECTION 

2)  PROBABILITY  OF  DETECTION 

3)  PROBABILITY  OF  DETECTION 

4)  COMPARSION  OF  M-S  CURVES 

vs.   S/N 
vs.   Pfa 
vs.    N 

Select  a  menu  number: 

Figure  8:  Submenu  Screen 

From  the  submenu  the  user  will  have  another  set  of 
choices.  He  can  choose  the  item  he  wants  to  study.  After  he 
chooses  one  of  the  items,  the  following  selected  screen  will 
be  seen  in  Figure  9,  Figure  10,  Figure  11  and  Figure  12.  The 
user  can  follow  the  instructions  on  the  screen  to  key  in  the 
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arguments  he  wants  to  study.  Then  a  data  input  screen  seen  at 
Figure  13  will  display  the  data  for  double  check.  After  the 
completion  of  above  procedure,  a  selected  graphic  output  will 
appear  on  the  screen  as  in  Figure  14 .  After  the  graphics 
display,  the  user  can  press  'enter'  to  get  the  next  screen  as 
shown  in  Figure  15.  At  this  point,  the  user  can  either  choose 
to  go  back  to  the  main  menu  or  exit  to  print  out  the  graphic 
display.  A  selected  second  submenu  and  third  submenu  are  shown 
in  Figure  16  and  Figure  17;  users  can  follow  the  same 
procedure  to  choose  the  items  they  want  to  study. 


%********************************************************* 

%      THIS  PROGRAM  RETURNS  THE  PLOTS  FOR  THE  NUMBER 

%  OF  PULSES  AND  SWERLING  CASE  SPECIFIED  IN  THE  PARAMETERS. 

%  THE  PLOTS  WILL  BE  STORED  IN  METAFILES  UNDER  THE  NAME 

%  "RADAR. MET"  FOR  AN  EASY  PRINT  OUT. 

% 

%   (A)    the  swerling  case  number  has  to  be  determined  now 

%********************************************************* 

echo  off 
Enter  the  case  number  you  want  to  study 


Figure  9:  Subprogram  Descriptions  Screen  (a) 


*********************************************************** 
(B)  .    The  number  of  radar  pulses  the  program  is  to 

integrate  needs  to  be  an  integer  between  1  and  600 
*********************************************************** 

cho  off 

Number  of  Pulses  to  be  integrated  is  n  =  10 


Figure  10:  Subprogram  Descriptions  Screen  (b) 
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****************************************** ****************************** 
(C) .    The  probability  of  false  alarm  rate  curves  (pfa)  to  be  plotted 
must  now  be  determined.   Each  choice  of  a  pfa  will  result  in  a  different 
curve  being  plotted  on  the  graph.   You  need  to  choose  the  following; 

1.  The  smallest  pfa  curve  to  be  plotted,  pfamin  =  ? 

2.  The  largest  pfa  curve  to  be  plotted,  pfamax  =  ? 

3.  The  step  size  between  pfamin  aand  pfamax,  pfastep  =  ? 

If  you  wish  to  plot  only  one  curve  then  enter  the  same  value  for 
pfamin  and  pfamax. 

The  suggested  default  step  size  to  use  is  that  of  PFASTEP  =  10, 
which  is  quite  sufficient.   It  is  suggested  that  pfamin  and  pfamax 

be  powers  of  10  as  that  is  the  normal  choice. 

***********************************************************************; 


Figure  11:  Subprogram  Descriptions  Screen  (c) 


The  signal  to  noise  ratio  (S/N)  in  dB  for  which  you  wish  to 
plot  needs  to  be  determined.   The  choices  you  need 
to  make  are; 

1.  The  smallest  S/N  point  to  be  plotted,  sdbmin  »  ? 

2.  The  largest  S/N  point  to  be  plotted,  sdbmax  -  ? 
Remember  that  S/N  must  be  entered  in  dB. 

3.  The  steps ize  between  sdbmin  and  sdbmax  -  ? 

******************************,***********. *****  ** 


Figure  12:  Subprogram  Descriptions  Screen  (d) 
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\    CHECK  YOUR  PARAMETERS  ! 
% 


echo  off 

The  case  number  you  is     1.00 

The  number  of  pulses  you  choice  are     1.00 

The  max  false  alarm  probability  you  choice  is  1.00e-12 

The  min  false  alarm  probability  you  choice  arel.00e-12 

The  pfa  stepsize  is     10.00 

The  max  signal-to-noise  ratio  you  choice  is    10.00 

The  min  signal-to-noise  ratio   -10.00 

The  S/N  stepsize  is     1.00 

% 

%IF  THE  PARAMETERS  ARE  CORRECT  PRESS  1 

%IF  THE  PARAMETERS  ARE  NOT   CORRECT  PRESS  2 

% 


Figure  13:  Input  Argument  Checking  Screen 
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Figure  14:  Selected  Result 
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If 

you  want  to  go  to  the  same  submenu 
or 

ENTER  CHOICE  =  1 

If 

you  want  to  go  to  the  main  menu 
or 

ENTER  CHOICE  =  2 

TO 

exit  this  whole  program 

PRESS  RETURN 

Figure  15:  Selected  menu 


SUB   MENU  [  NUMERICAL  CALCULATION]  

1)  CALCULATION  OF  DETECTION  PROBABILITY 

2)  CALCULATION  OF  THRESHOLD  LEVEL 


Figure  16:  Selected  second  menu 


SUB  MENU  —  [THE  DETECTION  RANGE  ANALYSIS]  

1)  RANGE  DETECTION  CURVES  WITH  DIFFERENT  PROBABILITY  OF  FALSE  ALARM 

2)  RANGE  DETECTION  CURVES  COMPARSION 

3)  RANGE  DETECTION  CURVES  USE  THE  DETECABILITY  FACTOR 


Figure  17:  Selected  third  menu 
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B.   Algorithms  for  M-S  models 

The  function  programs  in  Figure  5  are  based  on  the  Marcum 
and  Swerling  detection  theory. 

1.   Function  program  algorithms 

The  function  programs  can  be  executed  independently  as 
a  normal  MATLAB  function.  The  input  arguments  and  output  data 
are  listed  in  Table  [1]  at  the  end  of  this  chapter.  All  these 
programs  were  implemented  with  the  help  feature  of  the  MATLAB 
environment.  Typing  "help  <  name  of  the  function  >"  will 
explain  how  to  use  them  independently. 
a.  Threshold  computation 

The  detection  programs  start  with  the  calculation  of 
threshold  yb   by  applying  the  equation  (30)  in  Chapter  2. 

Pf  =PfAN,Yh)  =1-I(^,N-1)  =Y^Le~Yb  (85) 


fit  t=o   Kl 


The  false-alarm  probability  can  be  represented  as  a 
function  of  N  and  yb.  To  compute  probability  of  detection  for 
all  the  five  cases  Yb   is  required  for  a  given  pfa  . 

Since  Eq(85)  is  a  finite  power  series,  in  order  to  find  Yb 
for  given  Pfa,  a  recursive  computation  method  has  to  be  used 
rewriting  Eq(85) 
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Pfa[N,  Yb)  -g  lle-Y»=Pfa(N-l,  Yb)  t^e1'** 
=Pfa(N-l,Yb)  +L(N-1) 


(86) 


=pfa(#,  sg  +l(jv)  <87> 


=Pfa(tf,*g+L(tf-l)^ 


From  Eq  (86)  and  (87) 


L{n)  =L(N-1)^  (88) 

N 

For  a  single  hit  (N=l)     the  false  alarm  probability  and  the 
relation  with  the  following  term  (N=2)    are 

Pfa(l,Yb)=er»        ,        L(l)=Ybe-Y»  (89) 

Therefore  this  relation  allow  each  term  of  the  expansion 
to  be  based  on  the  value  of  the  proceeding  term,  therefore  an 
algorithm  can  be  formed  to  compute  the  values  of  the  detection 
threshold  yb   and  the  number  of  integrated  pulses  N. 

The  first  procedure  employed  in  the  algorithm  is  to  define 
an  empirical  threshold  level  Y0t  then  compute  P(N,y0)  .  On  the 
basis  of  this  empirically  determined  value  y0,  an  empirical 
suggested  Yb   is  given  by 
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Yb=N-sfN+2.3jZ(JL+J7J-l)  ,       where    L=-logPfa      (90) 

This  value  is  used  as  the  starting  point  to  compute 
P(N,y0)  then  compare  it  against  the  desired  value  of  Pfa  and 
the  difference  between  YN  and  YN+1.  This  value  of  correction  Ay 
can  be  calculated  by  using  Newton-Raphson  method  by  noting  hat 

ln_pfa(tf,3g 


e  Y"YJ!  V(iV-l)  ! 

The  procedure  is  repeated  until  the  correction  magnitude 
Ay/  (Ay+y)  is  indicated  that  Yb  is  within  a  sufficient 
accuracy.  A  computer  independent  algorithm  notation  is  given 
on  Figure  18.  In  Figure  18,  the  arrow  implies  a  specification. 
The  normal  execution  of  statements  is  carried  out  line  by 
line,  starting  at  the  top,  but  a  branch  may  be  designated  by 
an  arrow  which  results  from  the  execution  of  a  statement.  A 
conditional  branch  is  denoted  by  a  colon  statement,  and  the 
branch  is  executed  if  the  comparison  condition  specified  on 
the  arrow  is  satisfied.  Otherwise,  the  next  statement  in  the 
sequence  is  executed.  Notice  in  figure  18  that  the  program  is 
terminated  when  the  value  Ay/  (Ay+y)  is  less  or  equal  to  £ . 
The  value  of  £  can  be  assigned  to  be  10"6  or  10"12.  This 
accuracy  should  be  sufficient  for  application. 
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Figure  18:  Algorithmic  Program  for  Detection 
Thresholds/  yb 


The  MATLAB  function  for  executing  this  process  is  named 
prob.m  and  thresh. m  respectively.  In  Figure  18  e  is  the 
smallest  acceptable  tolerance  value.  Ay/ (Ay+y)  is  compared 
with  e,  and  when  it  is  less  than  or  equal  to  e,  the 
computation  will  stop  and  return  the  value  of  false  alarm 
probability.  In  thresh. m  the  smallest  tolerance  value  was  set 
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to  be  10"6,  this  value  is  sufficient  for  desired  accuracy. 

This  recursive  method  will  be  used  in  all  Swerling  cases 
to  find  the  threshold  level  yb_ 

b.    Swerling  case   1   algorithm 

Following  equation  is  used  to  compute  probability  of 
detection  for  case  1. 


PD  =1-J[     ~h    ,  (N-2)  ]  +(1+JL)       j( £ ,N-2)  e  ^m 

v/2FT  Wx  [n-(i/tfx)v/3FT] 

=Pfa  (tf-1,  Yb)  +  (1+-4:)  *"*  [l-^/a  W"1'  — T-  )1  e"^7^" 

iVx  1+-==. 

#x 

(92) 

It  is  obvious  to  see  that  for  N=l  the  detection  probability  is 
e-(Yb  /  i*nx>  _  rp^g  aig0rithm  in  computer  independent  notation  which 
uses  equation  (92)  to  compute  PD1  is  shown  on  Figure  19.  Notice 
in  Figure  18  that  the  detection  threshold  (Yb)  is  independent 
of  the  target  fluctuation  characteristics  so  that  the 
algorithm  given  in  Figure  18  is  used  to  determine  Yb  for  all 
target  types.  A  MATLAB  source  code  to  compute  probability  of 
detection  of  case  1  targets  is  given  in  Appendix  B.  The  name 
of  this  file  is  swerll.m. 
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Enter   Yb,N,X 


XN  i 

Y  i 
N 

A  < 

YM  i 

YMS  < 

M  < 

M  < 

N   -  1 

YM  < 

YMS  < 

Pz  * 
1 

A  f 

Pi  <■ 

Y  f 

Ax  ♦. 

A2  <■ 

Pj  <" 

P  <- 

-P  *• 


N    *    X 

Yb 
1 
0 

EXP    (-Y) 
YM 
0 

M   +    1 
M 

YM   o    Y/M 
YMS    +   YM 
YMS 
A 
1 

P2 
.  Y/(l    +1/XN) 

(1  +i/xII)N-1 

EXP    (-Yb/(1    +   XN)  ) 
Al    *   A2    o     (1    -    P2) 


Pi     *    Pi     

EXP  [  -Yb/(Xn  +    1)] 


•>   Exit  pDl 
>   Exit   Pol 


Figure   19:   Algorithmic  Program  for  Swerling  1  Target ,PD1 


c.    Swerling  case  2  algorithm 
Eq(63)    can  be  modified  as: 
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PD=1-I( 2 — 


=Pf*  (N, 


,N-1) 


(93) 


1+Nx 


A  computer  independent  algorithm  to  implement  the  above 
equation  is  shown  in  Figure  20. 


Enter   Yb 

,N 

Y 

<- 

Yb 

X 

<- 

(1    +   N    *    X) 

Xo 

<r- 

Y/X 

YM 

<- 

exp(-X0) 

YMS 

<— 

YM 

M 

«- 

0 

— ►m 

«- 

M    +    1 

ZZ.«.  M 

M 

X0*YM/M 
YMS    +    YM 
XMS       

YM 

<r- 

Xi'lO 

^-P«, 

i — 

->    Exit    P 

d2 


Figure  20:  Algorithmic  Program  for  Swerling  2  Target,  PDa 


50 


d.  Swerling  case   3   algorithm 

From  Eq(75) ,  the  equation  for  probability  of  detection  for 
case  3  is  as  follows: 


pD3=e     2     (1+^)2       folN=1 


Y%-2e-Y»    2Yb  Nx+2    N.2    ^±  (94) 


=  *-+Pfa(N-l,Yb)+(mZ±)       e** 
(N-2)   !  Nx+2  D  NX 

Nx  Nx+2  Nx+2 


An  algorithm  for  above  equation  is  presented  in  computer 
independent  notation  in  Figure  21. 
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Enter   Yb 

N,X 

XN 

C3 

N#X 

Y 

*— 

Y„ 

N 

1 

A 

*— 
<— 

0 

EXP    (-Y) 

YM 

YMS 

*— 

YM 

M 

f— 

0 

M 

<— 

M   +    1       -* 

N-l 

«— 

M 

YM 

<— 

YM      Y/M 
YMS    +    YM  ' 

YMS 

«— 

■P2 

♦- 

YMS 

1 

A 

A 

«— 

1 

PI 
Y 
"Al 

«— 

P2 

(1    *    2/X*)"'1 

A2 

<— 

EXP       (-2Yb/(X*+2)  ) 

A3 

<~ 

1-2    (N-2) / (X„   f 

2Yb/(XM   ♦    2) 

P3 

«— 

Al    *    A2    *    A3    *    (1 

-    P2) 

B 

1 

2 

N 

L 

*— 

N-2 

C 

«— 

L 

L 

♦— 

L    -    1 

0 

: 

L 

C 

♦— 

(C*L) 

B 

*— 

C 

P4 

^~ 

2ybH-i    evb/B.(Xw+2) 

p 

f_ 

PI      ♦     P9     +     P?        ___ 

_    i 

G 

*— 

EXP       -2Yb/(XN+2) 

H 

*- 

Z    XN   Yb/ (XN+2)2 

P 

*- 

G    (1    ♦    H)     

>Exit      P„, 

Figure  21:  Algorithmic  program  for  Swerling  3  Target,  P 
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e.    Swerling  case  4   algorithm 
Swerling    case    4    computer   program  modified   expression   can  be 
obtained  by   applying   the  power   series    in  Eq    (30)    into   Eq(81) 


as 


p„-i-  .    *!      ..T.   2 


n     4    H-PW+k,  l+x/2  )1 


D«      ~  (l+x/2)»h  *■  (N-k)  ! 

(-4= 


y    ta 


2H  i*iSE  ~7~        Sf  2 


•i-(-A^»"[Ee  -^ -£  km/wo- 

1  + 


Nx^       fa  fro 

2 


+ 


«  x*—    2        y^  Z        i 

2-   e  Si  £$   iC!  (JT-K)  ! 

y 

i+M        fa  ml        fa  K\(N-K)\ 

2 

( —  )m 

+  (!+**)"£  e     '-   I—] 

2  ^  ml  J 


m-2N 


( —  )* 

2AT-1  *  ZEf'         1+- 


-E 


e 


i*^=  2 


(95) 


( —  )IB 

-(      y    )    \^Nx  Nx 

-V   e           2     2          ry^           Nl  ,         2        \K/         1         \ 

jfo                            Ail           l£*  JC!  (tf-JO  !  l+Nx/2'       W2' 


JT      xa 


-(      r    )      14.AX 
i*iS  ~ 


"P-(W'TTW2»\Ee        '   — H 


m-N 


[1-Y  £! ( ? )JC, 1 )N.,n 

fa   K\  {N-K)  !      l+Nx/2         1+NX/2 
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with  this  expression  the  detection  probability  computation  for 
this  case  therefore  can  be  simplified  according  to  equation 
(95)  and  avoid  the  computation  of  infinite  power  series.  A 
computer  independent  algorithm  is  at  Figure  22. 


Enter  Yb 

N,X 

xN 

<r- 

X 

M 

<- 

0 

MK 

<— 

(2/ (XN  +2) )N 

Y 

<- 

Yb 

B 

«- 

2Y/ (XN  +  2) 

YM 

<— 

e  -B 

YMS 

YM 

M  +  1 

M 

N 

YM 

<— 

YM  o  B/M 

YMt* 

4 

YMSS  +  YM 

— - —  SUM 

<— 

YMS 

MKS 

<— 

MK 

—  YM 

<— 

YM  o  B/M 

SUM 

<— 

SUM  +  YM  (1-MKS) 

MK 

<— 

XN  o  MK  o  (2N-M)/2(M-N+1) 

MKS 

f- 

MKS  +  MK 

M 

<- 

M  +  1 

-£.2N 

: 

M 

P4 

«- 

SUM    >Exit   Pd4 

Figure  22:  Algorithmic  program  for  Swerling  4  Target/  pM 


54 


f.    Swerling  case   5  algorithm 

In  case  5,  the  expression  of  Marcum  steady  target  model  in 
Eq  (38)  can  also  be  modified  by  substituting  the  infinite 
power  series  for  the  modified  Bessel  function  IN(X)     given  by 

(-)2K 
IWU)  =(-)*£  KKN+X)  ! 

and  interchanging  the  order  of  summation  and  integration,  one 
obtain  the  probability  detection  for  case  5  as: 

Pos'e-"±  ^  Tf  e-^]  (97) 

tf=0    A>      m-0       /"* 


By  interchanging  the  order  of  the  summation  results  in  a 
efficient  computation  representation  as: 

N-l  vm       °°       v/n     m-N  .        .  K 

P=y  e-Y>^-+Y  e'Y^  (l-y  e-^Mf)       (98) 

D5  h      ">*■  fa     ™!     k        ^ 

The  right  hand  side  of  equation  (98)  has  two  terms,  the  latter 
term  is  a  infinite  summation  of  power  series,  therefore  a 
recursive  evaluation   is  necessary  for  computing  PD5. 

The  method  to  handle  this  infinite  power  series  is  to 
separate  above  equation  into  two  terms .  Let 

where  L  represent  the  integrated  pulses  where  L  represent  the 
integrated  pulses .  Let  L   z  N,    then  PL   can  be  represented  as 
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Af-1       -,rm       L  rrin  m-N 


L  k  *\    fa         m\         fa  K\ 

the  error  term  P„  can  be  represented  as 

P  =  f  e-r-g(l-?e-*(Jpf)Jr)  (lOD 

if  let  PL  =PD5  then  the  truncated  error  term  will  be  Pe 


(102) 


Let  L  be  large  enough  such  that 

«L     yK  L+2~N    (NX\k~-Nx  (103) 

m=0     A"       Jc^O       *"' 

Therefore 

PD5-PL*e  (104) 

An  upper  bound  can  be  found  to  limit  the  desired  accuracy 
and  avoid  the  unnecessary  computation.  The  computer 
independent  algorithm  is  in  Figure  23 . 
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enter   Yb/N,S 


XN 

<— 

X  *  N 

YM 

«— 

EXP  (-Yb) 

YMS 

<h- 

YM 

<— 

1 

YM 

<— 

YM  °  Yb/M 

YMS 
M 

<— 
<- 

YMS  +  YM 

rl  +  X 

—  SUM 

<— 

YMS 

XB 

f- 

EXP  (-XN) 

XBS 

<— 

XB 

— «■"  YM 

<— 

YM  o  Yb/M 

YMS 

f- 

YMS  +  YM 

SUM 

<— 

SUM  +  YM  (1  -  XBS) 

M 

<— 

M  +  1 

K 

<r- 

M  -  N 

XB 

<- 

XB  o    XN/K 

XBS 

<— 

XBS  +  XB 

<€ 

. 

(1  -  YMS) (1  -  XBS) 
SUM   

P 

<- 

>    eXIT  WITH  pD5 


Figure  23:  Algorithmic  program  for  Swerling  5  Target,  P 
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2.   Solving  for  round  off  error 

Round  off  error  can  take  place  during  the  calculation  of 
case  5  .  That  is  because  during  the  calculation  of  probability 
of  detection,  if  the  signal-to-noise  ratio  x  is  zero,  the 
output  detection  probability  will  become  the  probability  of 
false  alarm,  then  the  criterion  in  Eq  (103)  will  become 

L^N    (Nx)ke'Nx  (105) 

s(i-i)  (l-  T    {m.te     ) 

k-0  K- 

£0 

Since  we  can  not  set  the  minimum  acceptable  error  to  be 
zero,  the  computation  process  will  not  terminate.  Due  to  this 
reason,  trade-off  between  the  accuracy  and  the  maximum  value 
of  output  data  should  be  made.  An  algorithm  to  compute  the 
threshold  yb  to  deal  with  the  round  off  error  is  in  Figure  24. 
The  corresponding  MATLAB  file  name  is  Marcum.m  which  can  be 
used  to  compute  the  maximum  acceptable  threshold  level  yb  when 
the  value  of  signal-to-noise  ratio  is  zero  and  the  number  of 
pulses   is  large. 
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p  -* 

0 

0 

LOG(Y) 

EXP(-Y) 

M+l 

N-l 

LOG(M) 

YLX+YLOG-SUNK 

TSUM+EXP(YLX-Y) 

Y+LOG  (XPrA*TSUM)  *TSUM/EXP  ( YLX-Y) 

e 

Yl >  Yb 

*i     ■ 

0 
0 
LOG(Y) 

P 

M+l 

LOG(N) 

YLX+YLOG-SUNK 

EXP (YLX-Y) 

Yl    >   Yb 


Figure  24:  Algorithmic  Program  for  Detection  Thresholds,  Yb 
with  round  off  error  preventation 
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3.   Solving  for  underflow  and  overflow 

As  mentioned  before,  the  maximum  acceptable  number  for 
MATLAB  is  nearly  e-709-782712893384040999999,  the  underflow  and  overflow 
problem  arise  due  to  the  fact  that  the  power  series  in  M-S 
equations  has  the  form 

e-rXl        ;   e-i^  (106) 

m\  ml 

If  the  calculation  for  this  power  series  is  over  the 
maximum  or  below  the  minimum  number  that  MATLAB  can  represent, 
then  the  MATLAB  will  return  a  string  named  NaN  or  empty  matrix 
( [  ] )  and  the  whole  computation  process  will  be  meaningless . 
If  the  computation  needed  large  integrated  pulse  number  and 
large  threshold  or  signal-to-noise  ratio,  a  special  effort 
has  to  be  made  to  retain  the  accuracy.  One  method  is  let  the 
series  be  of  the  form 

k  k 

e-y  l-=exp[-y+kln(y)  -J^ln(N)]  =e"P        (107) 

and  let  the  value  e"p  be  compared  to  the  smallest  acceptable 
value  e-709-782712893384040999999  in  MATLAB,  of  course  e-709-782712893384040999999 
is  much  smaller  than  the  tolerance  number  €  assumed  in  the 
program.  If  pis  greater  than  709.7821...  then  increase  the 
number  K  in  above  equation  until  it  is  less  than  709.7821... 
and  then  start  the  summation  process.  This  method  will  save 
MATLAB  computation  of  power  series  with  a  large  number  and 
avoid  the  under  flow  problems.  An  algorithm  for  using  this 
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method  to  compute  the  detection  probability  is  in  Figure  (25)  , 
and  the  MATLAB  source  code  is  throshm.xn  and  is  given  in 
Appendix  B . 
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Figure  25:  Algorithmic  porgram  for  computing  the  probability 
of  dtection,PdS  with  underflow  preventation 
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4.   Input  arguments  and  output  data 

The  table  below  summarizes  all  input  and  output  data  for 
the  programs  . 


PROGRAM 

INPUT   DATA 

OUTPUT    DATA 

THRESH. M 

Pfa  N 

Yb 

PROB .  M 

Yb  ,N 

Pfa 

SWERL1.M 

Pfa'N,        X 

Pdi 

SWERL2 . M 

Pfa/N,       X 

PD2 

SWERL3 . M 

Pfa/N#       X 

PD3 

SWERL4.M 

Pfa/N,        X 

Pd4 

SWERL5 . M 

Pfa/N,        X 

PD5 

THRESHM . M . m 

N,    P£a 

Yb 

MARCUM.M 

Yb    ,    N   ,     x 

Pdi       or    Pfa 

POINT.  M 

N,     x 

Pd    (numerical) 

RADPOINT.M 

N,     x 

Pd    (output   data) 

RADAR. M 

Pfa/N,        X 

CHART   S/N  vs    Pd 

RADAR1.M 

Pfa/N,        X 

CHART   Pfa  vs    Pd 

MEYER . M 

Pfa/N,        X 

CHART  N  vs    Pd 

RANGE. 1 

Pfa/N,        X       , 

HART   R  vs    Pd 
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RANGE. 3 


Pfa,N,   x  , 


CHART  R  vs  Pd 


Table  1  j  Input  arguments  and  output  data 


Where 

N=Number  of  pulses  to  be  integrated 
X  =average  signal-to  -noise  ratio 

Pd=  probability  of  detection 
pfa=probability  of  false  alarm 
R=  the  desired  detection  range 
ns=Swerling  case  number 
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IV. RESULTS 

The  M-S  detection  probability  curves  can  be  obtain  by 
applying  the  MATLAB  programs  in  Chapter  3  .  The  results  can  let 
the  user  evaluate  the  detection  probability  detection 
performance  easily.  The  detection  range  curves  can  illustrate 
relationship  between  the  detection  probability  and  the 
detection  range  by  input  the  properly  detection  loss. 
A.   PROBABILITY  OF  DETECTION  CURVES 

As  mentioned  previously,  the  M-S  model  arise  from  the 
different  fluctuations  of  target  cross  section.  MATLAB 
programs  can  be  used  to  plot  probability  of  detection  versus 
per  pulse  signal-to-noise  ratio  for  given  probability  of  false 
alarm  for  five  cases .  It  can  also  be  asked  to  determine  per 
pulse  signal-to-noise  ratio  required  for  given  PD  and  Pfa. 
This  required  signal-to-noise  ratio  can  be  used  to  compute 
maximum  detection  range  from  the  radar  range  equation. 

Fig  26  and  Fig  27  compare  the  five  target  models  for  a 
false  alarm  of  10~6,  and  the  number  of  integrated  pulse  as  N=10 
and  N=100  respectively.  When  the  probability  is  larger  than 
0.33,  all  four  cases  in  which  the  target  cross  section  is  not 
constant  requires  greater  signal-to-noise  ratio  than  the 
constant  cross  section  case.  This  increase  in  signal-to-noise- 
ratio  will  cause  a  reduction  in  detection  range.  Therefore  if 
the  characteristic  of  the  target  cross  section  are  not 
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properly  take  into  account,  the  actual  performance  of  the 
radar  might  not  measure  up  to  the  performance  which  is 
predicted  from  the  constant  target  cross  section.  A  greater 
signal-to-noise  ratio  is  required  when  the  fluctuations  are 
correlated  pulse  to  pulse  (case  1  and  case  3  )  than  when  the 
fluctuations  are  uncorrelated  pulse  to  pulse (case  2  and  case 
4)  .  Fig  27  also  indicates  that  when  the  number  of  integrated 
pulses  is  large,  the  case  2  and  case  4  will  approach  to  the 
constant  target   case  (case  5) . 
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Figure  26:  Probability  of  detection  curves  for  five  target 
models 
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Figure  27:  Probability  of  detection  curves  for  five  target 
models 

Since  both  the  false  alarm  rate  Pfa  and  detection 
probability  PD  are  specified  by  the  system  requirement,  the 
radar  designer  computes  required  signal-to-noise  ratio  from  M- 
S  curves  and  uses  this  to  compute  the  maximum  detection  range. 
The  greater  the  number  of  pulses  integrated,  the  greater  is 
resulting  overall  signal-to-noise  ratio.  This  results  in 
greater  probability  of  detection,  but  it  will  require  longer 
scan  time.  Figure  28  and  Figure  29  give  plots  of  Pd  VS  SNR 
for  selected  Pfa's  for  case  1  for  N  equal  to  100  and  500 
respectively. 
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Figure  28:  Probability  of  detection  curves  for  Swerling  case 
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Figure  29:  Probability  of  detection  curves  for  Swerling  case 
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Figure  30  and  Figure  31  are  the  Swerling  case  2  plots  of 
Pd  Vs  SNR  with  10  and  100  pulses  integrated  respectively.  From 
Figure  30,  signal-to-noise  ratio  of  4.8  db  approximately  is 
required  to  yield  a  probability  of  detection  of  0.8  with  10 
pulses  integrated  and  probability  of  false  alarm  of  10"6. 


Figure  30:  Probability  of  detection  curves  for  Swerling  case 
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Figure  32  and  Figure  33  show  the  Swerling  case  3  with  10 
and  100  pulses  integrated  respectively. 
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Figure  32:  Probability  of  detection  curves  for  Swerling  case 
3 


Figure  33:  Probability  of  detection  curves  for  Swerling  case 
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Figure  34  and  Figure  35  show  Swerling  case  4  pulse  to 
pulse  fluctuation  with  10  and  100  pulses  integrated 
respectively. 
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Figure  34:    Probability  of  detection  curves  for  Swerling  case 
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Figure  35:  Probability  of  detection  curves  for  Swerling  case 
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Figure  3  6  and  Figure  37  are  Swerling  case  5  steady-state 
target  with  10  and  100  pulses  integrated  respectively. 
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Figure  36:   Probability  of  detection  curves  for  Swerling  case 
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Figure  37:   Probability  of  detection  curves  for  Swerling  case 
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B.   THE  DETECTIOON  RANGE  CURVES. 

The  ASR-9  radar  is  used  as  a  sample  to  compute  maximum 
detection  range.  ASR-9  radar  is  designed  for  surveillance  of 
the  terminal  area  around  airports.  A  significant  feature  of 
this  S-band  radar  is  its  use  of  moving  target  detection  (MTD) 
processor  which  provides  coherent  integration  among  other 
features.  The  parameters  of  the  ASR-9  radar  are  given  in  table 
below.  Also  tabulated  is  the  empirically  determined  search 
detection  range  on  a  1  m2  target  for  the  various  Swerling 
models  with  a  PD  =0.9  and  Pfa  =10"6. 


ASR-9  PARAMETERS  AND  RANGE  CALCULATION 

Pt-12oo  kw        ranges 

£-2900  MHz  Rx  -39.41  nm 

T-1.05  us  R,  -52.41  nm 

a  -1  Rj-46.68   nm 

Gt-33.5  db  R4-57.29   nm 

G£-33.5  dB  R,-62.63   nm 
NF-5  dB 
L-12dB 

6-1.3°  Losses 

9-75°/sec 

PRF-1200  pps  Transmitter   2  dB 

B^pp-150  Hz  Receiver      2   dB 

Po-0.9  Mismatch      1  dB 

Pfa-10'*  Integrator    1  dB 

p-l  Collapsing    1  dB 

Beam  Shape    3  dB 
MTI  .  2   dB 


Total         12  dB 


Table  2  :   ASR-9   PARAMETERS  AND  RANGE  CALCULATION 
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Figure  38:  Probability  of  detection  curves  for  five  target 
models 


Figure  39:  Range  detection  curves 
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C.   Numerical   data  output   from  probability  of   detection 
programs . 

A  table  of  values  below  is  included  for  the  purpose  of 
checking  the  programs . 


N 

Y„ 

cast 

x=3.  162278 

x*10 

x*31. 62278 

x»100 

x-316.2278 

1 

13.8155055 

1 

.0361810926687 

.2848035870497 

.65475592794524 

.8721557722127 

.957383959769 
4 

2 

.0361810926687 

.2848035870497 

6547559279452 

.8721557722127 

.957383959769 

4 

3 

.0202712698585 

2918843939447 

.77944622707425 

.9625961846032 

.995941497060 
5 

4 

.0202659173728 

.2918820913596 

.77944704886892 

.9652566475676 

.995941410232 

7 

5 

.0045852017166 

.2480487019045 

.99722505883065 

.9999995192496 

.999999519246 
9 

3 

19.1291681 

1 

.1971745188135 

.5769062596383 

.83647025963833 

.9446918645347 

.982126063254 
8 

2 

.1630815180577 

.7468903012970 

.97820931643215 

.99901694621S7 

.999965068946 

1 

3 

.1784365463374 

.6869964903258 

94510656525289 

.9933271492948 

.999288753543 
2 

4 

.1368031960489 

.8340491394154 

.9977036152469 

.9999940489495 

.999999991764 
8 

5 

.0888130400938 

.9727252372625 

.99999920934189 

.9999992093418 

.999999209341 
8 

N 

V. 

can 

x»3. 162278 

x»10 

x-31. 62278 

x«100 

x-316.2278 

10 

32.7013405 

1 

.  .4855434894507 

.7911151202538 

.9282416527967 

.9765960615308 

.992532970604 

2 

2 

.73389703346785 

9989667533193 

99999988582414 

.9999999999973 

1 

3 

.56937472704103 

.9139452495156 

.9800436879910 

.9988084719731 

.999877746289 

5 

4 

.78178938715194 

.9999629142286 

.99999999999768 

1 

1 

5 

.8355167894356 

.9999991524271 

.99999915242839 

.9999991524283 

.999999152428 

3 

N 

Y„ 

cut 

x-3. 162278 

x«10 

x«31  62278 

xalOO 

x»316.2278 

30 

63.5481801 

1 

.69852634183006 

.8917070566484 

.96429073448000 

.9885538115323 

.99639654382453 

2 

.99944729942719 

.9999999999989 

1 

1 

1 

3 

.8147880711125 

.9759019998584 

99728588213879 

.99971801636709 

.99997145717056 

4 

.999931272763883 

1 

1 

1 

1 

5 

.999998744645914 

.9999992846098 

.99999928460987 

.99999928460987 

.99999284460987 
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V.   CONCLUSIONS  AND  RECOMMENDATIONS 

A.  CONCLUSIONS 

Radar  detection  theory  is  outlined  in  this  thesis  and 
computer  programs  are  developed  for  probability  of  detection 
calculations  in  MATLAB .  These  programs  are  based  on  radar 
detection  theory  as  developed  by  Marcum  and  Swerling.  These 
programs  give  more  accurate  results  than  the  commonly 
available  programs  based  on  the  detectability  method  which  is 
empirical  is  native.  A  user's  guide  and  instruction  are 
included  in  the  thesis.  Simple  examples  of  the  programs  usage 
are  illustrated  in  the  thesis. 

Students  can  use  the  programs  written  for  this  thesis  to 
investigate  radar  system  performance.  These  programs  are  cost 
effective,  convenient  to  use  and  easy  to  reproduce  since  they 
are  run  on  personal  computers  that  are  readily  available  to 
students.  In  particular,  the  program  can  aid  students  by 
removing  a  major  computational  burden  to  allow  him  to  perform 
real  world  detection  calculations.  These  programs  can  also  be 
used  for  calculations  in  radar  development  work. 

B .  RECOMMENDATIONS 

Due  to  an  overflow  underflow  problem  the  programs  are 


76 


limited  to  integration  of  600  pulses  in  the  probability  of 
detection  calculations  for  all  5  cases.  Most  of  the  time  this 
does  not  present  any  limitation  as  most  radars  integrate  fewer 
than  600  pulses.  However  it  may  be  desirable  to  remove  this 
limitation  by  employing  a  Chernoff  bound  or  other  methods. 
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Appendix  A 

The  gram-Charlier  Series  is  a  series  expansion  of  a 
probability  density  function  in  terms  of  the  Gaussian 
function,  its  derivatives  and  the  moments  of  the  original 
density  function. 

If  P(x)  is  a  probability  density  function  that  is  nearly 
Gaussian  and  x  has  zero  mean  and  unit  variance,  then  the  Gram- 
Charlier  series  expansion  of  P(x)    is  given  by 


P(x)  =V]  a£  <J)(i)  (x)  where 

(A-l) 

<|)(0)  (x)=4>(x)=-J^e-*2/2 

v/2¥  (A-l. a) 


A(i)(x)-dftM 

dx1 


Coefficient  a£     can  be  evaluated  in  terms  of  the  Hermite 
polynomials : 


<l)i(x)=_UJLle(-x2/2)H,  (x)  (A-2) 

Some  typical  polynomial  are 

tf0(x)=l,  ff4(x)  =x4=6x2+3 

/^  (x)  =x,  H5  (x)  =x5-10x3+15x 

tf2  (x)  =x2-l,  H6  (x)  =x6-15x4+4  5x2-15 

tf3  (x)  =x3-x,  Hn  (x)  =x7-21x5+105x3-105x 

A  recursive   relation   for    these  polynomials    is 


(A-3) 
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Hi+1  (x)  =xHi  (x)  -iH^  {x)  (A-4 ) 

A  particularly  useful  fact  is  that  the  hermite  polynomials  and 
the  derivatives  of  the  Gaussian  function  are  biorthogonal . 
That  is, 

f'°y(x)Hj(x)dx 

=  (-Di  i!  6i:,  (A-5) 

.((-DM!,    i=j 

0        i*j 

To  evaluate  coefficient  ait  both  side  of  Eq(A-l)  are 
multiplied  by  H£(x)  and  result  is  integrated  from  -°°  to  °°. 
With  Eq(A-3),the  result  simplifies  to 

■1) 


ai=  [°°p(x)Hi(x)    dx  (A-6) 

1   \  J  -oo 


The  first  three  coefficients  a0    a1    and  a2    can  be  found  by 
substituting  Eq(A-3)into  Eq(A-6) 


aQ  =  (~P(x)dx=l 

a^-f"  xp(x)dx=0  (A-7) 

J      —00 


With  these  results,  the  Gram-Charlier  series  expansion  of  p  (x) 

in  Eq(A-l)  can  be  simplified  to 

when  the  mean  of  random  variable  x    is  not  zero  and/or  its 
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p(x)=^(x)+J^^(x)  (A-8) 


1^3 

variance  o  is  not  unity,  the  Gram-Charlier  expansion  of  p(x) 
is  obtained  by  expanding  a  related  density  function  g(t)  in  a 
Gram-Charlier  series,  where 

t=^l£  (A-9) 

a 

By  a  simple  transformation  of  variables  it  follows  that 


dx 
dt 


a         a 


(A-10) 


If  the  Gram-charlier  series  expansion  for  git)    is 

co 


i-0 


then  from  Eq(A-lO)  and  Eq(A-ll)  ,  the  expression  for  p  (x) 
becomes 


pW'-E^t— )  (A"12) 


If  both  sides  of  Eq(A-12)  are  multiplied  by  H7[(x-  x    )    I    o] 
and  the  result  integrated  from  -«  to  »,  coefficient  cL  can  be 
evaluated  with  relation  (A-6)  as 
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c  ±^lrp{x)H ^ (**>«*  (A-i3) 

1  !      J  —  O 

Evaluating    Eq(A-13)     for    the    first    few   coefficients    ct   yields 


c0  =  l,  c4  =  — (a4-3) 


1_ 

4! 


c^c^Q,  c5  =  --^rya5-l0a3  (A-14) 


3Ta'  C6  =  l6~!) 


c3  =  -^a,  cfi  =  ^-  (afi-i5a4+30) 


where  o^  is  the  ith  central  moment  of  p  (x)  ,    normalized  by  ct : 

a^—T   (x-x)ip(x)dx  (A-14) 

The  04  can  be  expressed  in  terms  of  the  conventional  moments 
of  the  distribution  p(x);  thus  if  m„  denotes  the  nth  moment 
of  p (x)    that  is 

mn=  [~xnp(x)dx  (A-15) 

J   —  « 

then  06,  through  c^  are  related  to  noncentral  moments  n^   by 

/n3-3/772/r?1+2/n1 

"^ ^ 

2  2  (A-16) 

_  .m5  - 5^^ + 1  Qm^ml  - 1  On^m?  +4/??! 

a5  - 

a5 

_  m6-6m5m1+l 5m^ni1  -20m2m1  +15m2m1  -5/% 
6  a6 

In    summary,     to    obtain   a   Gram-Charlier    series    expansion, 
the     moments     m„     are     found     from    Eq(A-16)  ,      or     they     can    be 
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conveniently  computed  from  the  characteristic  equation  of  the 
distribution.  Next  the  oq  are  evaluated  using  Eq(A-17)  .  The  oq 
are  then  substituted  into  Eq(A-14)  to  obtain  coefficient  cL  of 
Gram-Charlier  expansion. 
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APPENDIX  B 

%FILE  NAME:  MS . M 

clear 

clg 

clc 

echo  on 

%THIS  PROGRAM  IS  DESIGNED  TO  ALLOW  THE  STUDENT  TO 

%VARY  THE  PARAMETERS  OF  THE  VARIOUS  SWERLING  CASES 

%IN  ORDER  TO  STUDY  THE  EFFECTS. 

%  CASE  #:     DESCRIPTION 

%  1.  Returned  pulses  are  of  a  constant  amplitude  over  one 
scan,  but  are  uncorrelated  from  scan  to  scan. 

%  2 .  Returned  pulses  are  uncorrelated  from  pulse  to  pulse 
and  correlated  from  scan  to  scan. 

%  3 .  Returned  pulses  are  of  a  constant  amplitude  over  one 
scan,  but  are  uncorrelated  from  scan  to  scan. 

%  4 .  Returned  pulses  are  uncorrelated  from  pulse  to  pulse 
and  correlated  from  scan  to  scan. 

%  5.  The  static  case  with  constant  S/N  and  pulse  amplitude 

echo  off 
pause 

clc 

k=menu( 'MAIN  MENU' , 'THE  M-S  CURVES 

' NUMERICAL  DETECTION  PROBABILITY  CALCULATION 

'RANGE  DETECTION  CURVES') 

if  k«l 

ms curve 

elseif  k==2 

number 

elseif  k==3 

rmenu 

elseif  k==4 

density 

end 

%   Now  restart  the  process  with  some  different  choices 

%  Clear  the  workspace  of  unnecessary  data  to  avoid 

conflicts 

clear  i;  clear  sav;  clear  sdb;  clear  w,  clear  pfa,  clear 

p; clear  ns; 

clear  i;  clear  savl;  clear  sdbl;  clear  w,  clear  pfal,  clear 

P; 

%  Check  and  see  if  any  info  is  to  change 

clc 

echo  on 

%   If  you  want  to  restart  the  process       ENTER  CHOICE  =  1 

% 

%       or 
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% 

%   To  exit  this  programPRESS  RETURN 

echo  off 

choice  =  input ( 'CHOICE  =  ') 

if  choice==l 

ms 
else 
end 
clc 

echo  on 

%   You  have  chosen  to  exit  the  program 
% 

%    Bye! 
echo  off 


%  FILE  NAME:  RADAR. M 

clc 

format  long 

echo 

%*  +  +  +  +  +  +  +  +  +  +  +  +  +  +  *  +  *  +  +  +  +  +  +  +  *  +  +  +  ■*■*•*  +  ■*•*  +  ■*■*•*■*■*  +  ■*■*■*  +  +  +  ■*■*■*■■*■*•*•*•*■*•* 

%      THIS  PROGRAM  RETURNS  THE  PLOTS  FOR  THE  NUMBER 

%  OF  PULSES  AND  SWERLING  CASE  SPECIFIED  IN  THE  PARAMETERS. 

%  THE  PLOTS  WILL  BE  STORED  IN  METAFILES  UNDER  THE  NAME 

%  "RADAR. MET"  FOR  AN  EASY  PRINT  OUT. 

% 

%   (A)    the  swerling  case  number  has  to  be  determined  now 

echo  off 

ns=  input ('Enter  the  case  number  you  want  to  study  '); 

echo  on 

%  *******•*****■**■*■**  +  ■*•*****•****  +  *-*•  +  +  •*■*  +  *  *********************** 
* 

%   (B) .    The  number  of  radar  pulses  the  program  is  to 
%         integrate  needs  to  be  an  integer  between  1  and 
600. 

echo  off 

n=input ( 'Number  of  Pulses  to  be  integrated  is  n  =  '); 

clc 

echo  on 

%   (C) .    The  probability  of  false  alarm  rate  curves  (pfa) 

to  be  plotted 

%   must  now  be  determined.   Each  choice  of  a  pfa  will  result 

in  a  different 

%   curve  being  plotted  on  the  graph.   You  need  to  choose  the 

following; 


84 


%    1.   The  smallest  pfa  curve  to  be  plotted,  pfamin  =  ? 

% 

%    2.   The  largest  pfa  curve  to  be  plotted,  pfamax  =  ? 

% 

%    3  .   The  step  size  between  pfamin  aand  pfamax,  pfastep  = 

-> 

% 

%    If  you  wish  to  plot  only  one  curve  then  enter  the  same 

value  for 

%       pfamin  and  pfamax. 

% 

%       The  suggested  default  step  size  to  use  is  that  of 

PFASTEP  =  10, 

%       which  is  quite  sufficient.   It  is  suggested  that 

pfamin  and  pfamax 

%       be  powers  of  10  as  that  is  the  normal  choice. 

echo  off 

pfamin  = input ( 'pfamin  =  ') 
pfamax  =input ( 'pfamax  =  ')• 
pf astep=input ( 'pfastep  =  ' ) 
clc 
echo  on 

%   (D) .    The  signal  to  noise  ratio  (S/N)  in  dB  for  which 

you  wish  to 

%         plot  needs  to  be  determined.   The  choices  you  need 

%  to.  make  are ; 

% 

%    1.   The  smallest  S/N  point  to  be  plotted,  sdbmin  =  ? 

% 

%    2.   The  largest  S/N  point  to  be  plotted,  sdbmax  =  ? 

% 

%  Remember  that  S/N  must  be  entered  in  dB. 

% 

%       3 .      The  steps ize  between  sdbmin  and  sdbmax  =  ? 

echo  off 

sdbmin  =input('In  dB  the  sdbmin  is   ') 

sdbmax  = input ('In  db  the  sdbmax  is   ') 

sdbstep=  input ('In  db  the  sdbstep  is  ') 

clg 

clc 

echo  on 

% 

%  CHECK  YOUR  PARAMETERS  ! 

echo  off 

fprintf('The  case  number  you  is  %8.2f \n' ,ns) ; 

fprintf('The  number  of  pulses  you  choice  are  %8.2f\n',n); 
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fprintf ( 'The  max  false  alarm  probability  you  choice  is 

%8.2e\n' ,pfamax) ; 

fprintf ('The  min  false  alarm  probability  you  choice 

are%8 . 2e\n ' ,pfamin) ; 

fprintf ('The  pfa  stepsize  is  %8 .2f \n ' , pf astep ' ) ; 

fprintf ('The  max  signal-to-noise  ratio  you  choice  is 

%8 . 2f \n' , sdbmax) ; 

fprintf ('The  min  signal-to-noise  ratio  %8 . 2f \n ' , sdbmin) ; 

fprintf ('The  S/N  stepsize  is  %8 .2f \n ' , sdbstep ' ) ; 

echo  on 

% 

%IF  THE  PARAMETERS  ARE  CORRECT  PRESS  1 

%IF  THE  PARAMETERS  ARE  NOT   CORRECT  PRESS  2 

% 

echo  off 

choice=input ( 'CHOICE=  '  )  ; 

if  choice==2 

radar 

elseif  choice==l 

end 

echo  on 

%   Now  run  the  calculations  for  the  above  data  and  plot  it 

echo  off 

%delete  radar. met 

format  long 

nmax  =  700; 

axis  (  '  square  '•)  ; 

w=[  sdbmin  sdbmax  0   100]; 

axis  ( w)  ; 

if  n  <  1  |  n  >  600  |  n-=fix(n), 

error (' Number  of  pulses  input  exceeds  dimensions') 
end 
if  ns  <  1  |  ns  >  5  |  ns~=fix(ns), 

error (' Swerling  case  input  does  not  correspond  to 
allowable  choices ' ) 
end 

pf a=pfamax; 

sdb  =  sdbmin : sdbstep: sdbmax ; 

kk=length(sdb) ; 

sav  =  10  . "(sdb/10) ; 

if  ns  ==  1 

%       swerlingl 
while  pfa  >=  pfamin 
for  i=  l:kk 

%i  =  1 : sdbmax- sdbmin+1 
p(i) =swerll(pfa,n, sav(i) ) ; 
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i 


p(i)=100*p(i) ; 
end 

plot(sdb,p, ' - ' ) ; grid; title ( ' SWERLING  1  '  )  ; 
xlabel ( ' (S/N) 1,  signal-to-noise  ratio,  DB ' ) , 
ylabel (' Probability  of  detection  in  percent  ') 
text(sdbmin  +  l/20*sdbmax, 95  ,  ['n  =  ' ,  int2str (n)  ] ) 
text(sdbmin  +  l/20*sdbmax,  90, ['pfa  = 
1 , num2str (pf amax) ,  '  to  ' , num2str (pf amin) ] ) 

text(sdbmin  +  l/20*sdbmax,  85, ['pfastep  = 
1 , num2str (pfastep) ] ) ; 

gtext ( [num2str (pfa) ] ) ; 
hold  on 

pfa  =  pfa/pfastep; 
end 
pause 

meta  radar 
elseif  ns==2 
%swerling2 
pf a=pfamax; 
while  pfa  >=  pfamin 
for  i=  l:kk 

%  for  i  =  1  :sdbmax-sdbmin+l 
p(i) =swerl2 (pfa,n, sav(i) ) ; 
p(i)=100*p(i) ; 
end 

plot(sdb,p, '-' ) , title ( 'SWERLING  2' ) ; 
grid,  xlabel (' (S/N) 1,  signal-to-noise  ratio,  DB ' ) , 
ylabel (' Probability  of  detection  in  percent') 
text(sdbmin  +  l/20*sdbmax, 95  ,['n  =  ' , int2str (n) ] ) 
text(sdbmin  +  l/20*sdbmax, 90, [ 'pfa  =  ' ,num2str (pfamax) , ' 
to  ' , num2str (pfamin) ] ) 

text(sdbmin  +  l/20*sdbmax, 85,  [ 'pfastep  = 
' , num2str (pfastep) ] ) 

gtext ( [num2str (pfa) ] ) ; 
hold  on 

pfa  =  pfa/pfastep; 
end 
pause; 
meta  radar 
elseif  ns==3 
%swerling3 
if  n  =s  1 

error ( 'algorithm  for  swerling  3  does  not  work  for  n=l ' ) 
end 

while  pfa  >=  pfamin 
for  i=l:kk 

%for  i  =  1  :sdbmax-sdbmin+l 
p(i) =swerl3 (pfa,n, sav( i) ) ; 
p(i)=100*p(i) ; 
end 
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plot  (scab,  p,  '  -  '  )  ,  title  (  'SWERLING  3  '  )  ; 
grid,  xlabel ( ' (S/N) 1,  signal-to-noise  ratio,  DB ' ) , 
ylabel (' Probability  of  detection  in  percent') 
text(sdbmin  +  l/20*sdbmax, 95  ,['n  =  ' , int2str (n) ] ) 
text(sdbmin  +  l/20*sdbmax,  90  ,  [  'pf  a  =  '  ,  num2str  (pf  arnax)  , 
to  ' , num2str (pfamin) ] ) 

text(sdbmin  +  l/20*sdbmax, 85 ,[ 'pf astep  = 
' , num2str(pf astep) ] ) 

gtext ( [num2str (pf a) ] ) ; 
hold  on 

pfa  =  pfa/pfastep; 
end 
pause; 
met a  radar 
elseif  ns==4 
%swerling4 

while  pfa  >=  pfamin 
for  i=l:kk 

%for  i  =  1 :sdbmax-sdbmin+l 
p( i) =swerl4 (pfa,n, sav(i) ) ; 
p(i)=100*p(i) ; 
end 

plot(sdb,p, ' - ' ) ; title ( 'SWERLING  4  '  )  ; 
grid;  xlabel (' (S/N) 1,  signal-to-noise  ratio,  DB '  )  ; 
ylabel ( 'Probability  of  detection  in  percent') 
text(sdbmin  +  l/20*sdbmax, 95  ,['n  =  ' , int2str (n) ] ) 
text(sdbmin  +  l/20*sdbmax, 90, [ 'pfa  =  ' ,num2str (pfamax) , 
to  ' , num2str (pfamin) ] ) 

text(sdbmin  +  l/20*sdbmax, 85, [ 'pf astep  = 
' ,num2str(pf astep) ] ) 

gtext ( [num2str (pfa) ] ) ; 
hold  on 

pfa  =  pfa/pfastep; 
end 
pause 

met a  radar 
else 
%swerling5 

while  pfa  >=  pfamin 
for  i=l:kk 

%  for  i  =  1 :sdbmax-sdbmin+l 
p(i) =swerl5 (pfa,n, sav(i) ) ; 
p{i)=100*p(i) ; 
end 

plot(sdb,p, '-' ) , title ( 'SWERLING  5' ) ; 
grid,  xlabel (' (S/N) 1,  signal-to-noise  ratio,  DB ' ) , 
ylabel (' Probability  of  detection  in  percent') 
text(sdbmin  +  l/20*sdbmax, 95  ,['n  =  ' , int2str (n) ] ) 
text(sdbmin  +  l/20*sdbmax, 90, [ 'pfa  =  ', num2str (pfamax) , ' 
to  ' ,num2str (pfamin) ] ) 

text(sdbmin  +  l/20*sdbmax, 85 , [ 'pfastep  = 
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1 , num2str (pfastep) ] ) 

gtext ( [num2str (pfa) ] ) ; 
hold  on 

pfa  =  pf a/pf astep; 
end 
pause 

meta  radar 
end 
hold  off 


%FILE  NAME:  SWERL1.M 

function  pd=swerll (pfa, n, sbar ) 

%  RETURNS  THE  VALUE  OF  THE  PROBABILITY  OF  DETECTION  FOR  THE 

SWERLING1 

%  CASE,  GIVEN  THE  NUMBER  OF  PULSES  n,  THE  PROBABILITY  OF 

FALSE  ALARM 

%  PFA  AND  THE  AVERAGE  SIGNAL  TO  NOISE  RATIO  S/N. 

%  EXAMPLE:  SWERL1 ( PFA, N, SBAR) ; 

format  long 

yb=thresh(pfa,n) ; 

y=yb/ (l+n*sbar) ; 

cte= (1+1/n/sbar) ; 

pd=prob(n-l,yb)+cteA (n-1) *exp(-y) * (l-prob(n-l/yb/cte) ) ; 

return 

end 


%  FILE  NAME:  SWERL2 .M 
% 

%  RETURNS  THE  VALUE  OF  THE  PROBABILITY  OF  DETECTION  FOR  THE 

SWERLING2 

%  CASE,  GIVEN  THE  NUMBER  OF  PULSES  n,  THE  PROBABILITY  OF 

FALSE  ALARM 

%  PFA  AND  THE  AVERAGE  SIGNAL  TO  NOISE  RATIO  S/N. 

%  EXAMPLE :  SWERL2 ( PFA , N , SBAR ) ; 

■^  ************************************************************ 

function  pd=swerl2 (pfa, n, sbar) 

yb=thresh(pfa,n) ; 

y=yb/ (1+sbar) ; 

pd=prob(n,y) ; 

return 


%  FILE  NAME:  SWERL3 .M 

%  RETURNS  THE  VALUE  OF  THE  PROBABILITY  OF  DETECTION  FOR  THE 

SWERLING3 

%  CASE,  GIVEN  THE  NUMBER  OF  PULSES  n,  THE  PROBABILITY  OF 
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FALSE  ALARM 

%  PFA  AND  THE  AVERAGE  SIGNAL  TO  NOISE  RATIO  S/N. 

%  EXAMPLE:  SWERL3 ( PFA, N, SBAR) ; 

%** ********* ******* ***** 

function  pd=swerl3 (pf a, n, sbar ) 

yb=thresh (pf a, n) ; 

y=l; 

for  k=l:l:n-2; 

y=yb/k*y; 

end; 

yl=2*yb/ (n*sbar+2) ; 

y=y*exp(-yb) *yl; 

y2=( (n*sbar+2) /n/sbar) A (n-2) *exp(-yl) * (1-2* (n-2) /n/sbar+yl) ; 

y2=y2* (l-prob(n-l,yb*n*sbar/ (n*sbar+2) ) ) ; 

pd=y+prob ( n- 1 , yb ) +y2 ; 

return 

%  FILE  NAME:  SWERL4.M 
^***********»*************w*********»*********************** 

%  RETURNS  THE  VALUE  OF  THE  PROBABILITY  OF  DETECTION  FOR  THE 

SWERLING4 

%  CASE,  GIVEN  THE  NUMBER  OF  PULSES  n,  THE  PROBABILITY  OF 

FALSE  ALARM 

%  PFA  AND  THE  AVERAGE  SIGNAL  TO  NOISE  RATIO  S/N. 

%  EXAMPLE:  SWERL4 ( PFA, N, SBAR) ; 
*•*•*•*••*•*****•**••***•*••***•*•••*•***•**•**••****•••**»• 

function  pd=swerl4 (pfa, n, sbar ) 

yb=thresh(pfa,n) ; 

v=2*yb/  (sbar-t-2)  ; 

ZK=(2/ (sbar+2) ) An; 

YM=exp(-v) ; 

YMS=YM; 

for  M=l:n-1 

YM=YM*v/M; 

YMS=YMS  +  YM; 
end 

SUM=YMS ; 
ZKS=ZK; 
for  M=n: (2*n-l) 

YM=YM*v/M; 

SUM  =  SUM  +  YMM1-ZKS); 

ZK  =  sbar*ZK*(2*n-M) / (2*(M-n+l) ) ; 

ZKS  =  ZKS  +  ZK; 
end 

pd=SUM; 
return; 

%FILE  NAME : SWERL5 . M 

%  RETURNS  THE  VALUE  OF  THE  PROBABILITY  OF  DETECTION  FOR  THE 
SWERLING5 
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%  CASE,  GIVEN  THE  NUMBER  OF  PULSES  n,  THE  PROBABILITY  OF 

FALSE  ALARM 

%  PFA  AND  THE  AVERAGE  SIGNAL  TO  NOISE  RATIO  S/N. 

%  EXAMPLE :  SWERL5 ( PFA , N , SBAR ) ; 

&  +  +  +  +  +  +  +  +  +  +  +  +  +  +  +  +  +  +  +  *  +  +  +  +  *  +  +  +  +  +  +  +  ■*  +  ■*■*•*  +  ■*•*  +  *■*  +  ■*■**■+*+■  +  •*■*  +  +  +  +  ■*  + 

function  pd=swerl5 (pf a, n, sbar ) 

yb=thresh (pf a, n) ; 

z=n*sbar ; 

YM=exp(-yb) ; 

OK  =  0; 

k=n; 

YMS=YM; 

for  i=l:k-l 

YM=YM*yb/i; 

YMS  =  YMS  +  YM; 
end 

SUM=YMS ; 
XB=exp(-z) ; 
XBS=XB; 
while  OK  ==  0 

YM=YM*yb/k; 

YMS  =  YMS  +  YM; 

SUM=  SUM+  YMM1-XBS); 

el=l-YMS; 

XB  =XB*z/ (k-n+1) ; 

XBS  =  XBS  +  XB; 

e2=l-XBS; 

er=el*e2; 

if  er  <=  le-6 
OK=l; 

end 

k  =  k  +  1; 
end 

pd  =  SUM; 
return; 

%  FILE  NAME:  THRESH. M 

%*  RETURNS  THE  THRESHOLD  yb  FOR  A  GIVEN  PROBABILITY  OF  FALSE 

ALARM  (pfa) 
%*  AND  A  GIVEN  NUMBER  OF  PULSES 
%*  EXAMPLE:  THRESH ( Pfa, N) . 

function  y=thresh(pfa, n) 

format  long 

if  n  <  1  |  n-=fix(n) , 

error ('Wrong  input  number  of  pulses'); 

break 
end 
if  n==l 

y=-log(pfa)  ; 
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return 
end 

1  =  -loglO(pfa) ; 
n2=sqrt(n) ;12=sqrt(l) ; 
y=n-n2+2.3*12* (12+n2-l) ; 
comp=le-6 ; 
ratio=l; 

while  comp  <=  ratio 
p=prob (n,  y) ; 
ym=l; 

for  k=l:l:n-l; 

ym=ym/k*y; 
end; 
ym=ym*exp ( -y ) ; 

dely= (p/ym)  *log(p/pfa) ;    %  correction  magitude 
y=y+dely; 

ratio=abs (dely/y) ; 
end 
return 

%  FILE  NAME:  PROB . M 

%*   THIS  PROGRAM  RETURNS  THE  FALSE  ALARM  PROBABILITY  P(n,yb) 
as  in       * 

%*   SUMMATORY(ybAk/fact(k) )  for  k=0  to  n-1,  AND  IS  A  COMMON 
FUNCTION     * 

%*   FOR  THE  PROBABILITY  OF  DETECTION  COMPUTATIONS 

* 

%*   EXAMPLE:  PROB(N,YB) 

• 

function  p=prob(x,y) 
if  x<  0  |  x~=fix(x), 

error (' Number  of  pulses  should  be  integer  and  greater  than 
zero ' ) , break, end 
inf ==1.797 693 13 4862069* 10~3 08 ; 
E=709. 7827 12893 3 84040999999; 
P=0; 
if  x==0 

return; 

elseif  x==l 

p  =  exp(-y) ; 

return; 

else 

%  •  [if  x>l] 

t=l; 
for  k=l:l:x-l; 
t=t/k*y; 
P=t+p; 
end 
p=(p+l) *exp(-y) ; 
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end 

return; 

%File  name:  THRESHM.M 

function  yb=threshm(pfa, n) 

tol=le-6; 

E=7 09. 7 827 12 8 93 3 840409999 99; 

1  =  -loglO(pfa) ; 

n2=sqrt (n) ; 12=sqrt (1) ; 

y=n-n2+2.3*12*(12+n2-l) ; 

plog=logl0 (pfa) ; 
xpfa=  1/pfa;' 
ok=0  ; 
while  ok==0 

if   y  <  E 

m=0; 
ylx=0; 

ylog=log(y) ; 
tsum=exp ( -y) ; 

for  m  =  0 : 1 :n-l 
m=m-»-l; 

sunk=log (m) ; 
ylx=  ylx+ylog  -sunk; 
tsum=tsum+exp (ylx-y) ; 
end 
yl=y+log  (xpfa*tsum)  *tsum/exp  (ylx-y)  ; 

if  abs(y-yl)  <  tol 

yb=y 1 ; 

ok=l; 

return 
else 

y=y  1  ; 

ok=0  ; 
end 


else 


m=0; 
ylx=0 ; 

ylog=log(y) ; 
okk=0 ; 
while  okk==0 


m=m+ 1 ; 

sunk=log (m) ; 
ylx=ylx+ylog-sunk ; 
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if  y-ylx  <=  E 

tsum=exp (ylx-y) ; 
okk=l; 
else 

okk=0; 
end 
end 
if  m  >=  n-1 

yl=y+log(xpfa*tsum)  *  tsum/exp( -y+ylx) ; 
if  abs(yl-y)  <  tol 
yb=y 1 ; 
ok=l; 
break; 
return 
else 

y=yl  ; 
ok=0  ; 
end 


else 


for  m=m:l:n-l 

sunk=log (m) ; 

ylx=ylx+ylog-sunk; 

tsum=tsum+exp( -y+ylx) ; 
end 

yl=y+log (xpf a*tsum) *tsum/exp ( -y+ylx) ; 
if   abs(yl-y)  <=  tol 

yb=yl ; 

ok=l; 

break ; 

return 
else 

y=yl  ; 

ok=0  ; 
end 


end 
end 
end 


%FILE  NAME:  MARCUM.M 

function  [p,ymo] =bound(n, y, x) 

format  long 

err=le-6; 

E=709. 7827128933 84040999999; 


%******    start  the  program 


•*•*••••••••*••••*•■*••••*••*•■•'*"*•**■*'**•* 


lamdal=(n/ (2*y) ) ~2+(n*x) /y; 

lamda=l- (n/ (2*y) ) -sqrt ( lamdal) ; 

auxl= (-lamda*y) + (n*x*lamda) / (1-lamda) ; 


94 


aux2=exp(auxl) / (1-lamda) ^n; 
if  aux2  <=  err 
aux3=n* (x+1) ; 
if  y  >=  aux3 

p=0 ; return 
else 

p=l ; return 
end 
else 

ib=n-l; 

k=0; 

x=n*x; 

if  x  <=  E 

xk=exp(-x) ; 

xks=xk; 

m=0; 

if  y<=  E 

ym=exp(-y) ; 

yms=ym; 
% „, 


for  m=l:n-l 

ym=ym*y/m; 
yms=yms-»-ym; 

if  err>(l-yms) 
p=yms; 
ymo=ym 
x=x/n; 

return 
end 
end 
m=n-l; 
sum=yms ; 
ymo=ym 
ok=0; 

while  ok==0 
k=k+l; 
m=m+l; 
ym=ym*y/m; 
yms=yms+ym; 
sum=sum+ym*  ( 1-xks )  ; 
xk=xk*x/k; 
xks=xks+xk; 
ans= ( 1-yms ) * ( 1-xks ) ; 
if  ans<=  err 
ok=l; 
end 
end 
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p = s  um ;  ymo  =ym  ; 
else 

ymly=0  ; 
ylog=log(y) ; 
okk=0; 

while  okk==0 
m=m+ 1 ; 
xmf =log (m) ; 
ymly=ymly+ylog-xmf ; 

if  (y-ymly)  <=  E 

ym=exp ( -y+ymly) ; 
okk=l; 
else 

okk=0; 
end 
end 

if   m  <=  ib 
yms=ym; 


for   ml=m:n-l 

ym=ym*y/m; 
yms=yms+ym; 
if  err>(l-yms) 

p=yms; 

ymo=ym 

x=x/n; 

ook=l; 

return 
end 

end 
sum=yms ; 
ymo=ym 
ok=0; 

while   ok==0; 
k=Jc+l; 
m=m+ 1 ; 
ym=ym*y/m; 
yms=yms+yin; 
sum=sum+ym*  ( 1-xrs )  ; 
xk=xk*x/k; 
xrs=xrs+xk; 
ans= (1-yms) * (1-xrs) ; 
if   ans   >    err 

ok=0  ; 
else 

p=sum; 
x=x/n; 
break; 
return; 


96 


ok=l; 
end 
end 

p = s  um ;  ymo =ym  ; 

else      % for  m>  ib; 

ymo=0 ; return 
yms=ym; 
sum=0 ; 
kd=m-n; 
if   kd  >  x 
P=0; 
x=x/n; 
break; 
return; 
else 

okkk=0; 

while  okkk==0; 
k«k*l; 
kk=xk*x/k; 
xks=xks+xk; 

while  (k-kd)  <  0 
okkk=0; 
end 
end 

ooook=0; 
while  ooook==0; 

sum=sum+ym*  ( 1-xrs )  ; 
xk=xk*x/k; 
xrs=xrs+xk; 
ans= ( 1-yms ) * ( 1-xrs ) 
while  ans  <  err 
p=sum; 
x=x/n; 
ooook=l; 
break; 
return 
end 
k=k+l; 
m=m+ 1 ; 
ym=ym*y/m; 
yms=yms+ym; 
ooook=0; 
end 

end 


end 
end 


else 

xklx=0 ; 
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xlog=log(x) ; 
ooook=0 ; 
while  oook==0 

k=k+l; 

xklx==xklx+xlog-xkf ; 

while  (x-xklx)  >  E 
oook=0 ; 

end 
end 

xk=exp ( -x+xklx) ; 
ib=ib+k; 

xks=xk; 

m=0; 

if  y<=  E 

ym=exp(  -y)  ; 
yms=ym; 


for  m=l:n-l 

ym=ym*y/m; 
yms=yms+ym;- 

if  err>(l-yms) 
p=yms; 
ymo=ym; 
x=x/n; 
return 
end 
end 
sum=yms  ; 
ymo=ym; 
ok=0  ; 

while  ok==0; 
k=k+l; 
m=m+ 1 ; 
ym=ym*y/m; 
yms=yms+ym; 
sum=sum+ym* ( 1-xrs ) ; 
xk=xk*x/k; 
xrs=xrs+xk; 
ans= ( 1-yms ) * ( 1-xrs ) ; 
if  ans<=  err 

ok=l; 
end 
end 
p=sum; 

else 

ymly=0  ; 
ylog=log(y) ; 
okk=0; 
while  okk==0 
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m=m+ 1 ; 

xmf =log (m) ; 

ymly=ymly+ylog-xmf ; 

if  (y-ymly)  <=  E 

ym=exp ( -y+ymly) ; 
okk=l; 
else 

okk=0; 
end 
end 

if   m  <=  ib 
yms  =ym  ; 


ook=0 ; 

while   ook==0 
while  m~=ib 
m=m+ 1 ; 
ym=ym*y/m; 
yms=yrns+ym; 
if    err>(l-yms) 
p=yms; 
ymo=ym; 
x=x/n; 
ook=l; 
break; 
return 
else 

ook=0; 
end 
end 
end 
sum=yms ; 
ymo  =ym  ; 
ok=0; 

while  ok==0; 
k=k+l; 
m=m+l; 
ym=ym*y/m; 
yms=yms+ym; 
sum=sum+ym* ( 1-xrs ) ; 
xk=xk*x/k; 
xrs=xrs+xk; 
ans= ( 1-yms) * ( 1-xrs) ; 
if  ans  >  err 

ok=0; 
else 

p=sum; 
x=x/n; 
break ; 
return; 
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ok=l; 

end 

end 


else      % for  m>  ib; 

ymo=0 ; 
yms=ym; 
sum=0  ; 
kd=m-n; 
if   kd>x 
p=0; 
x=x/n; 
break; 
return; 
else 

okkk=0; 

while  okkk==0; 
k=k+l; 
kk=xk*x/k; 
xks=xks+xk; 

while  (k-kd)  <  0 
okkk=0 ; 
end 
end 

ooook=0; 
while  ooook==0; 

s\im=sum+ym*  ( 1-xrs )  ; 
xk=xk*x/k; 
xrs=xrs+xk; 
ans= ( 1-yms ) * ( 1-xrs ) 
while  ans  <  err 
p=sum; 
x=x/n; 
ooook=l ; 
break; 
return 
end 
k=k+l; 
m=m+ 1 ; 
ym=ym*y/m; 
yms=yms+ym; 
ooook=0 ; 
end 

end 


end 
end 


end 
end 
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